{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### MCMC"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX0AAAD8CAYAAACb4nSYAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xt8VNW99/HPj3BVBBTQys0gcieY0EBQKmiFitUKWimXQ9XWCtSqx15oo/KIh4MeTulLbQu0cjg+4CMo1WoeTkGoCtqicidyswjYCEEfRBARuQhxPX9kZpwkk8yeZDJ7Lt/365UXM3vW3vObMPnttddaey1zziEiIpmhgd8BiIhI4ijpi4hkECV9EZEMoqQvIpJBlPRFRDKIkr6ISAZR0hcRySBK+iIiGURJX0QkgzT0O4DK2rRp47Kzs/0OQ0QkpWzcuPFj51zbaOWSLulnZ2ezYcMGv8MQEUkpZva+l3Jq3hERySBK+iIiGURJX0QkgyRdm75IJjp9+jSlpaWcPHnS71AkyTVt2pQOHTrQqFGjWu2vpC+SBEpLSznnnHPIzs7GzPwOR5KUc45Dhw5RWlpK586da3UMNe+IJIGTJ0/SunVrJXypkZnRunXrOl0RKumLJAklfPGirt8TJX0RkQyiNn2RoIkTv3r8xBP+xQEVY4kHD5/n4YcfZtGiRWRlZdGgQQOeeOIJCgoK4htHmCuvvJLf/OY35Ofn1/oY8+fPZ/LkyXTo0IFjx45x8cUXM3XqVC6//HIAHnzwQQYPHszQoUMj7l9UVES3bt3o1atXxNf/+Mc/ctZZZ3HLLbfEHO+RI0dYtGgRd955JwAffPAB99xzD88//3wtPmn8KOmLCG+99RZ/+ctf2LRpE02aNOHjjz/miy++8DssT0aPHs2sWbMAWLVqFTfddBOrVq2iZ8+eTJs2rcZ9i4qKuP766yMm/TNnzjBp0qRax3XkyBHmzJkTSvrt2rXzPeGDmndEBPjwww9p06YNTZo0AaBNmza0a9cOgGnTptG/f3/69OnDhAkTcM4B5TX1n/70pwwePJiePXuyfv16brrpJrp27cqUKVMAKCkpoUePHtx666307duXm2++mePHj1d5/7/+9a9cdtll9OvXj1GjRnHs2DEACgsL6dWrF3379uUXv/hF1M9x1VVXMWHCBObOnQvAbbfdFkq0lY/15ptvsmTJEiZPnkxubi579uzhyiuv5P7772fIkCH89re/5aGHHuI3v/lN6PhPP/00l19+OX369GHdunUAVcr06dOHkpISCgsL2bNnD7m5uUyePJmSkhL69OkDlHfc/+AHPyAnJ4e8vDxWrVoFlF+53HTTTQwfPpyuXbvyy1/+0ut/oWdK+iLCt771Lfbt20e3bt248847ef3110Ov3XXXXaxfv55t27Zx4sQJ/vKXv4Rea9y4MX/729+YNGkSI0aMYPbs2Wzbto358+dz6NAhAHbu3MmECRPYsmULLVq0YM6cORXe++OPP2b69Om88sorbNq0ifz8fB599FEOHz7Miy++yPbt29myZUvoRBJNv379+Mc//lFhW6RjXX755dxwww3MnDmT4uJiunTpApTX0F9//XV+/vOfVzn2559/zptvvsmcOXP44Q9/WGMcM2bMoEuXLhQXFzNz5swKr82ePRuArVu38swzz3DrrbeGRuQUFxezePFitm7dyuLFi9m3b5+nz+2Vkr6I0Lx5czZu3MjcuXNp27Yto0ePZv78+UB5k0lBQQE5OTmsXLmS7du3h/a74YYbAMjJyaF3795ceOGFNGnShIsvvjiUrDp27MigQYMAGD9+PKtXr67w3mvWrGHHjh0MGjSI3NxcFixYwPvvv0+LFi1o2rQpP/rRj3jhhRc466yzPH2W4JVIuFiONXr06GpfGzt2LACDBw/m6NGjHDlyxFNMla1evZrvf//7APTo0YOLLrqId999F4Crr76ali1b0rRpU3r16sX773uaR80ztemLCABZWVlceeWVXHnlleTk5LBgwQLGjBnDnXfeyYYNG+jYsSMPPfRQhTHiweagBg0ahB4Hn585cwaoOsSw8nPnHMOGDeOZZ56pEtO6det49dVXefbZZ5k1axYrV66M+jk2b95Mz549K2xr2LCh52OdffbZ1R470mdp2LAhX375ZWiblzH0kU5MQeG/x6ysrNDvMV5U0xcRdu7cya5du0LPi4uLueiii0IJrE2bNhw7dqxWHZF79+7lrbfeAuCZZ57hG9/4RoXXBw4cyBtvvMHu3bsBOH78OO+++y7Hjh3j008/5dvf/jaPP/44xcXFUd/r9ddfZ+7cudxxxx0Vtld3rHPOOYfPPvvM82dZvHgxUF5Tb9myJS1btiQ7O5tNmzYBsGnTJv75z39GPfbgwYNZuHAhAO+++y579+6le/funuOoC9X0RaLxYyhngoeMHjt2jLvvvpsjR47QsGFDLrnkEubOnUurVq244447yMnJITs7m/79+8d87J49e7JgwQImTpxI165d+fGPf1zh9bZt2zJ//nzGjh3LqVOnAJg+fTrnnHMOI0aM4OTJkzjneOyxxyIef/HixaxevZrjx4/TuXNn/vznP1ep6X/22WcRjzVmzBjuuOMOfve733k6oZ177rlcfvnlHD16lCeffBKA7373uzz11FPk5ubSv39/unXrBkDr1q0ZNGgQffr04dprr+UnP/lJ6Dh33nknkyZNIicnh4YNGzJ//vwKNfz6ZDVdZvghPz/faREV8UV1yT0BSf+dd96pkqjSQUlJCddffz3btm3zO5S0Eun7YmYbnXNRbyJQ846ISAZR0heRepOdna1afpJR0hcRySBK+iIiGURJX0Qkgyjpi4hkECV9EeHQoUPk5uaSm5vL1772Ndq3bx96Xl+zbW7atInly5dHfO2VV16hZcuW5OXl0a1bN4YMGcKyZctCr8+ePTt0c1MkK1euZM2aNdW+/uKLL4bmwxk/fjxFRUWe4/7yyy+ZMWNG6HlZWRlXXHGF5/39ppuzJLPFe976FNW6devQXaoPPfQQzZs39zSrZVBZWRlZWVkxveemTZvYtm0bw4cPj/j6VVddFUrGmzZt4sYbb+Spp55iyJAhFW50imTlypW0adOGgQMHVnntzJkz3HjjjTHFGi6Y9AsLC4HyqRL+/ve/1/p4iaaavkgKKtq8n0EzVtK5cCmDZqykaPP+enuv73znO3z961+nd+/ezJs3DyhPnK1atWLKlCkMGDCAdevWsWTJErp3784VV1zB3XffzciRI4Hyu31vu+02BgwYQF5eHv/zP//DiRMnmDZtGgsXLiQ3Nzfq3bD9+vXjgQceCM2bP2XKFB5//HEAHnvsMXr16sWll17K+PHj2bNnD/PmzWPmzJnk5uby5ptvMn78eH7+859z1VVXcf/99zNv3jzuvffe0PFXrFjBFVdcQbdu3XjppZcAqpQZPnw4q1evprCwkM8++4zc3FxuueWW0O8Cyk8IP/vZz+jTpw85OTmhz/XKK69w9dVXc9NNN9G9e3duueWWePzX1Iqnmr6ZDQd+C2QB85xzMyq9/jPgR8AZ4CDwQ+fc+4HXbgWCc6JOd84tiFPsIhmpaPN+7nthKydOlwGw/8gJ7nthKwAj89rH/f0WLFjAeeedx/Hjx8nPz+e73/0u55xzDp9++in9+vVj+vTpHD9+nG7duvHGG2/QqVMnvve974X2nzZtGsOHD2f+/Pl88sknFBQUsGXLFh588EG2bdsWSt7R9OvXj9///vdVtv/617/m/fffp3Hjxhw5coRWrVrxox/9iDZt2oSS9pw5c9izZw+vvvoqDRo0CJ28gvbt28frr7/Orl27GDp0aGgeoEhmzJjBvHnzQldG4ROiPffcc+zYsYO3336bgwcP0r9/fwYPHgyUX63s2LGD888/n4EDB7JmzZqIVyL1LWrSN7MsYDYwDCgF1pvZEufcjrBim4F859xxM/sx8GtgtJmdB0wF8gEHbAzs+0m8P4hIba0vOcyuA8cqbLu/cGno8SNr9wIwrqBTQuOqzswVO0MJP+jE6TJmrthZL0n/scceY8mSJQCUlpaGFgZp3LhxqJlkx44ddO/enYsuuggon4L4qaeeAsoXSHnppZdC7eAnT55k7969McdR3ZQxvXv3Zvz48YwYMSJ0dRHJqFGjaNAgcuPG9773PRo0aED37t3p2LFjhcnnYrF69WrGjRtHVlYWX/va1/jGN77Bhg0baNy4MQMHDuTCCy8EIDc3l5KSEl+SvpfmnQHAbufce865L4BngRHhBZxzq5xzweVw1gAdAo+vAV52zh0OJPqXgcgNeCIJll24lEVr91ZJ+NVZtHYv2YVLyQ47IfjhgyMnYtpeF6+88gp/+9vfWLNmDW+//TZ9+/YNzbzZrFmz0FTDNc3h5ZyjqKiI4uJiiouL2bt3b2hSslhEmjIZyptmJk2axLp168jPz6esrCzC3uk/ZbJXXpJ+eyB86ZbSwLbq3A68VMt9RepdPBK3n8m/XatmMW2vi08//ZTzzjuPZs2asX37dtavXx+xXO/evdm5cyf79u3DOReaghjgmmuu4Xe/+13o+ebNm4HYpjUuLi7mkUceqdKBW1ZWRmlpKd/85jeZOXMmBw8e5Pjx4zFPmfzcc8/hnOPdd99l3759dO3alezsbDZv3oxzjpKSEjZu3AiUz80PREzagwcP5tlnn6WsrIwDBw7wxhtv1Gnh9/rgJelbhG0RT2dmNp7yppzg2mCe9jWzCWa2wcw2HDx40ENIIrUT70TtR+KffE13mjWqOFKmWaMsJl8T//nYr7vuOo4fP86ll17KtGnTKCgoiFjurLPOYtasWQwdOpQrrriCdu3a0bJlSwCmTp3K8ePHQ6trPfTQQwB885vf5O233yYvLy9iR+6qVavIy8uje/fu3HPPPcyZM4chQ4ZUKHPmzBnGjRtH37596devH7/61a9CUzL/6U9/Ii8vjzfffDPq57zkkksYPHgw3/nOd5g7dy6NGzdmyJAhtG/fnpycHAoLC8nNzQ2Vv/322+nbt2+VDtmbb76ZHj16cOmllzJ06FAeffRRzj///Kjvn0hRp1Y2s8uAh5xz1wSe3wfgnPuPSuWGAr8HhjjnPgpsGwtc6ZybGHj+BPCac67qEjkBmlpZ6kukBP3I8lkRy94//K6IZarbPq74Jeoi1qmVizbvZ+aKnXxw5ATtWjVj8jXd66U9PxbHjh2jefPmOOeYOHEiOTk53H333b7GlK7qMrWyl9E764GuZtYZ2A+MAcZVerM84AlgeDDhB6wAHjGzcwPPvwXc5+E9ReKmaPN+7l0cfdWlusguXErJjOvq9T3Cjcxr73uSr+wPf/gDCxcu5NSpU+Tn51dZvUqSQ9Sk75w7Y2Z3UZ7As4AnnXPbzWwasME5t4Ty5pzmwHOBDpG9zrkbnHOHzezfKT9xAExzzh2ul08iEkGkhF9d7b6uEp34k83kyZOZPHmy32FIFJ7G6TvnlgHLKm17MOzx0Br2fRJ4srYBitRFbWv4LZpkcfRU5FEgNalL4nfOVRlFIlJZXVc71DQMkrZq08kaHIs/7t++Glm8KMYrg0W51341pt/j8opNmzbl0KFDtG7dWolfquWc49ChQzRt2rTWx1DSl9RWee6cQJKNNeHXdONV8LVFa73fULRo7d6Ybubq0KEDpaWlaPSaRNO0aVM6dOgQvWA1lPQl7cQz4Vcud38Mx120dm/FEQ81aNSoEZ07d47h6CK1ownXJK2EJ/xHls8K/URyQYsmMU+tEGt7fd+pkacOFvGLkr5kpAtaNOHqnhfUat9YEn9tOoNF6pOSvqSNWNrca5vwg2JJ/H7P1SMSTklf0kIsCT9es2XGchwlfkkW6siV1FOH1a5KZlwHE5fE7T3GFXRi0dq99XbDl0i8qaYvKc9rLb++7pb1WuNXbV+SgZK+pDSvCX/8wORYAOVf/ustv0OQDKekLxlh+sicej2+19r+G3s09ZT4S0lfUpbXWn6iljlUM4+kAiV9SWtdL2judwgiSUVJX1KS11p+/+zz6jmSilTbl2SnpC8pZ/G65GrWSZb3FfFCSV9STlndphNPGqrtix+U9CWlFDz8sqdyfq9gpdq+JCvdkSsp5cBnX0QtM6hLYtvxa1LdoupBmb7EoiSekr6kDK/NIQvvuKyeI/EmOEVDULQTgEgiKOlLWqm21lyH+XrqousFzdl14FiNZVTbl0RSm76khFTt9Ez0kFGRaFTTl7SRkM7TWlwxNMoyTkcZcqTaviSKavqS9DqnaC0/aFR+R79DEAlR0pek52VYfrIPkfQSX6o2YUlqUdKXpKapiEXiS0lfkpqXqYiTvZYfpNq+JAN15ErymjiRRwLj3DWuXSQ+lPQlaXmZSbNK7dmn8fheXdalNW/tOVRjGY3kkfqk5h2RBOrc5my/Q5AMp6QvSclL23aqtOVXZha9jNr2pb4o6Ysk2NgBqXmykvSgpC/JZeJEFuVeW2FyskhStZYfC9X2pT54SvpmNtzMdprZbjMrjPD6YDPbZGZnzOzmSq+VmVlx4GdJvAIXSWUlnyzhkeWzop7cROItatI3syxgNnAt0AsYa2a9KhXbC9wGLIpwiBPOudzAzw11jFfSnJcRO17axNPFsEdf8zsESTNeavoDgN3Oufecc18AzwIjwgs450qcc1uAL+shRpEK0qVN3EsT1a6PPk9AJJJJvCT99sC+sOelgW1eNTWzDWa2xsxGxhSdZBSvSyGKSO15SfqRLqZjWZq6k3MuHxgHPG5mXaq8gdmEwIlhw8GDB2M4tKQTL0shpttNS5qaQRLNS9IvBcLnhu0AfOD1DZxzHwT+fQ94DciLUGaucy7fOZfftm1br4cWEZEYeZmGYT3Q1cw6A/uBMZTX2qMys3OB4865U2bWBhgE/Lq2wUr68lKbTbdaflDJjOsifv7wkT3ZpO/nl8SKWtN3zp0B7gJWAO8Af3LObTezaWZ2A4CZ9TezUmAU8ISZbQ/s3hPYYGZvA6uAGc65HfXxQUREJDpPE64555YByyptezDs8XrKm30q7/cmkFPHGCXNeVkZS7VckfjQLJviu1hGBaSrYBNPTTdrafZNiQdNwyC+mlK01e8QRDKKkr746uk1tZgzP91MnAgTJ1LySfRZSjR8U+pKzTsiKaBCs4+aeKQOlPTFN8Faa03t2KFafpKviBUv4wo6RZ1/qMcDy/jHw99OUESSbtS8I5JiTpap61tqT0lffNF36vKoZdK+LV/EB0r64oujp8r8DiFpaT4eqU9K+iIiGURJXxIunRc9jxcvn19TUUttKOmLpCgvU1GLVKakL4mjRc9j0vWC5n6HIGlISV8kSfXPPi9qGXXoSqyU9CWpZNKi5yJ+UNKXhIl2pymkz6Ln8aLhmxJvmoZB4i98yoQnnvAvDhGpQklfEiK7cCmPhD2P1JmrDtzIvMzHM6VoK9NHar0iiU7NOyJpwMsU1SKgpC8J4GWhlMu6tE5AJCKipC/1zksttHObsxMQSeoKb/p6ZPms0E84deiKF0r6IiIZRElf6pXm2YkfL7+nYY++Vv+BSEpT0hdJI7s++tzvECTJKelLvXn1nQNRy6iWL5JYSvpSbw4cPeV3CGlHd+hKXenmLKlX0WbUFJHEUtKXeuFlnh017dROdXfohp9gLzHY/R/XJTIsSRFq3hFJQ2ec3xFIslLSl7h7cVNp1DKq5ddNs0b605Xa0TdH4u7E6S/9DiHt3divQ9Qy6tCVSJT0JeG0DKCIf5T0Ja681C69LAMo0Wn4ptSGkr6ISAbxlPTNbLiZ7TSz3WZWGOH1wWa2yczOmNnNlV671cx2BX5ujVfgknw0z07ilXyyJOKMm+H6Tl2ewIgk2UVN+maWBcwGrgV6AWPNrFelYnuB24BFlfY9D5gKFAADgKlmdm7dwxYRr46eKvM7BEkiXmr6A4Ddzrn3nHNfAM8CI8ILOOdKnHNbgMrDNq4BXnbOHXbOfQK8DAyPQ9ySZP7lv96KWka1/PpxQYsmfocgKcRL0m8P7At7XhrY5oWnfc1sgpltMLMNBw8e9HhoSSZv7DnsdwgZ6+qeF0Qtow5dCfKS9C3CNq/3+3na1zk31zmX75zLb9u2rcdDi4hIrLwk/VKgY9jzDsAHHo9fl30lRXjpKFTTTv3S8E3xykvSXw90NbPOZtYYGAMs8Xj8FcC3zOzcQAfutwLbJI2oo1AkdURN+s65M8BdlCfrd4A/Oee2m9k0M7sBwMz6m1kpMAp4wsy2B/Y9DPw75SeO9cC0wDZJExqmmTy83OlctHl/AiKRZOZpamXn3DJgWaVtD4Y9Xk95002kfZ8EnqxDjJICNG++//pnn8euA8dqLHPv4mJG5nkdhyHpSHfkSq0VPPxy1DKXdWmdgEgkKCvS0AmRMEr6UmsHPvsiapnObc5OQCQSNHqAOnSlZlo5S+qNap3JIbzp7f7hd/kYiSQDJX2pFS+1RS+1Tom/cQWdIEofS3bhUkpmaDnFTKTmHRGRDKKkLzHTMM3k52X4ptr2M5OSvkga0kI1Uh216UtMwodpVjc2v5F6cEWSlmr6EhMvwzRH5XeMWkbqn+bjkUiU9EVEMoiSvnimDtzUo9q+VKakL554majL1JQvkvSU9MWTexcXRy0zVjdjJSXV9iWckr6ISAbRkE2JSm35qe+yLq15a8+hKsNsw+fi6fHAMv7x8LcTHZokmGr6IhnAy2ynJ8u8Ln0tqUw1falRduHSqLM0qpafGszARcnrqu2nP9X0RTKEl4521fbTn5K+VMtLW36zRvoKpRIvU2T0nbo8AZGIX/QXK3VyY7+ISyNLkvIyRcbRU2UJiET8ojZ9iai6Wn54+36LZvr6pKJmjRpw4vSXNZYpePhl1j4wLEERSSKppi+1dn3fdn6HILXg5erMy8R6kpqU9KUKL235l3VpnYBIxE89HljmdwhSD5T0pVa8jPuW5OVlmK1G8qQnJX2pQLX8zOFlgrzwRXMkPagnTiKqblUsUC0/XYwd0In7opRR2376UU1fQrzU8rUUYuZR2356UdIXAIY9+pqncloKMb2UzLguahm17acXJX0BYNdHn0ctc0GLJgmIRBLt8dG5Uctovv30oaQvnjvrru55QT1HIn4Ymdfe7xAkgZT0xVNnnWbSTG9emnlU208PGr2T4TrrDzmzTZz41eNzb4ha/F/+6y0W3nFZPQYk9c1TTd/MhpvZTjPbbWaFEV5vYmaLA6+vNbPswPZsMzthZsWBnz/GN3ypKy9ddKrlZwYvbftv7DmcgEikPkVN+maWBcwGrgV6AWPNrFelYrcDnzjnLgEeA/4z7LU9zrncwM+kOMUtceDlct3LDTySHry27Wvq5dTmpaY/ANjtnHvPOfcF8CwwolKZEcCCwOPngavNlC6S2ZSirZ7KeVl4Q9KHl7Z9Tb2c2ry06bcH9oU9LwUKqivjnDtjZp8CwXv1O5vZZuAoMMU59/fKb2BmE4AJAJ06KckkwtNr9lbZVvkuXN2IJdXJLlzq6QQhycdL0o/0l1+5Kbi6Mh8CnZxzh8zs60CRmfV2zh2tUNC5ucBcgPz8fN0JUs8uuc9b561uxMpMJTOuY1HutaHnkdZFhvKrxekjcxIVlsSJl+adUiD8r78D8EF1ZcysIdASOOycO+WcOwTgnNsI7AG61TVoqYWJE0M/ZzycVrte0Lz+Y5Kk5WWBnEhXi5L8vCT99UBXM+tsZo2BMcCSSmWWALcGHt8MrHTOOTNrG+gIxswuBroC78UndKmNRWuj/6FmGfTPPi8B0Uiy8rpAjob8pp6oSd85dwa4C1gBvAP8yTm33cymmVlwYO9/A63NbDfwMyA4rHMwsMXM3qa8g3eSc05jvnzyzDpvNbPR6rwVvE2h7YCizfvrPxiJG083ZznnlgHLKm17MOzxSWBUhP3+DPy5jjFKnDgPzTpa91aCOrc5m7f2HIpa7t7FxZrKIYXoLzxDhDfr1DRXvta9lXDjCjpB2Peluk5djeZJHUr6GSC7cCmPeCinO28zXPiUDGGyDLzMrtx36nK2/NvwOAcl8aakL0D5H7ZIJKMHdApdKT5SQ61fN22lBs2ymea8zoyozlupiderQM3EmfxU009jXm/C0ph88aJRlnE6rJ2nulp/wcMvs/aBYQmNTbxTTT/dBG7AenXwSE83YYHG5Is3Xu/Q1mLqyU1JP00dOHrKUzl13kos1MyT+pT005CXu25BCV9qx+tEfEr8yUlt+ukgbKid1/H4WuRcamtUfkfPFQuN308+qumnEa9/iKBFzqVuKl8lPrJ8VuinMs3Pk1yU9NPEcxv2RS8UoGYdiQev3yNH+dq6khzUvJMG/vnx5xWG0tVECV/i6bIuravMzxNpKKfW1k0equmnuL5Tl3uaFEukPnRuc7bndZTVsZsclPSTXdjiJ5Vdct/SmG59Vy1f6kMs6yhnFy5VU4/PlPRTVN+pyz3ffAVK+FK/Yvl+vbHnsBK/j9Smn4IKHn5ZNXxJOuMKOlUZQVbdVA1q4/ePavopZtijr8V0m7sSviRSLN+37MKlWnXLB0r6KSS7cCm7Pvrcc3klfPFDLN+7excXq4M3wdS8kyIWrd0LYetT1HS3LSjhi7+iNfWEu3/4XbpzN4FU009yL24qjelOW1DCl+QQ6/dQzT2JoaSfxLILl3Li9Jcx7aOEL8kk1u/jvYuL6fHAsnqKRkBJP2nVpp1TCV+SUazfy5NlTu389cici2GwdwLk5+e7DRs2+B2GbwoefrnC6JxobfdBSviS7Lw0U1Zed1ft/N6Z2UbnXH60curITRKdC5dSm9Nvlml9W0kNkTp3ownW+LuefzYv/+zKeogq86h5x2fDHn2N7Fom/HEFnZTwJaWMK+jkeRGWcLs++pzswqVMKdpaD1FlFiV9nxRt3h/zuPtwas6RVDUqv2OtF/F5es1eLrl/mUb51IGadxJocd61BGdArtx2GQslfEl1wUV8Ym3uATjzpePexcXcu7iY9q2aMfma7ozMax/vENOWkn4CDHv0NXZ99DmPeGjDqanjtlGWMSq/YxwjE/HXuIJOvPrOAQ4cPVWr/fcfOcG9i4vZ8P5hpo/MiXN6l0r3AAAI9klEQVR06UlJvx4Ubd7PAy9u5fMvvE2KFm2Ejlls09eKpJJgrf+5DftqvGu3Jk+v2cvTa766arjgnMasfWBY/IJMI0r6cRKszXvldShmi2YNub5vu9qGJZIyRuV35J8ffx6XRYEOfPZFhbH+557ViKnf6a1mIDROP2ZFm/czc8VO9h85QZYZZVF+f16Te7jLurSmc5uzaxuiSFrw0uxTl76xsxtn8fCNOWlzItA4/TrqO3V51DnroyX8WKlWL/KVYLPPi5tKY56OpLp5/MN9/kVZqEMYyJhOYU9J38yGA78FsoB5zrkZlV5vAjwFfB04BIx2zpUEXrsPuB0oA+5xzq2IW/RhgjXwD46coF2E/7wpRVt5Zu0+ypwjy4yxBR0rdPyE7w/wsIcvTbja1OjDqXYvEtmN/TrUqbPXq/1HTnDfC1vZ8P5hVv3jIB8cOUGrsxpx8nRZhZPOuWc14rq+F4bKGFD5lBRsBUjGE0nUpG9mWcBsYBhQCqw3syXOuR1hxW4HPnHOXWJmY4D/BEabWS9gDNAbaAe8YmbdnHPel33yoGjzfu57YSsnTpcfNvifBzAyrz1TirZW6OQpcy70fPrInCr7J4LupBXxLljrh/IO39NlVa+yq6t4xdI5fOJ0GQvX7A3dLPnJ8dNVynxy/HSFfOKoemURbAWonIuqE63SGk9ebs4aAOx2zr3nnPsCeBYYUanMCGBB4PHzwNVmZoHtzzrnTjnn/gnsDhwvrmau2FklYZ84XcbMFTsBeGbtvoj7BbdH2r++tGjWUHfSitTBqPyOjCvoxLiCTpTMuI7HR+cS+z2+5Yk6+BMu3r2c4bkokmClc/+REzi+OlHU1w1oXpp32gPhWbMUKKiujHPujJl9CrQObF9Tad+4n76CTTJQ8Yz7QOBMXuZcxDa+4Nk4uH+sNQWvLmjRpEJNRUTiZOJERgIjAQI3LVZ3JeBVXf/eI+Wa8BxVxaSJ/K9A/2GwfPBEUR+1fS9JP9JJtPJvtLoyXvbFzCYAEwA6dYq9BtyuVTP2B36p4Zds7Vs1A8rb1yJdymWZVdg/1pEAwXa7ZG6/E8k0o8Ieh/fleWHUbURQddoFclEkP73qxxGvLmo8UdSBl6RfCoTfBtoB+KCaMqVm1hBoCRz2uC/OubnAXCgfsuk1+KDJ13Sv0ibfrFEWk6/pDsDYgo4V2uCCxhZ0rHb/6kTqBBaR5DR9ZE6Vv9Wizfu5/4UtHK80Iqh9q2Zc1aMtf964P67NveG5KJLwSmvl7fXBS9JfD3Q1s87Afso7ZsdVKrMEuBV4C7gZWOmcc2a2BFhkZo9S3pHbFVgXr+CDgjXr6jpCgv/p1Y3eiba/iKSPkXnta/zbzr/ovFAuSMTonWiV1njzdHOWmX0beJzyIZtPOuceNrNpwAbn3BIzawr8HyCP8hr+GOfce4F9HwB+CJwB7nXOvVTTeyX7zVkiIvEWj9E7Xm/O0h25IiJpwGvS13z6IiIZRElfRCSDKOmLiGQQJX0RkQyipC8ikkGU9EVEMoiSvohIBlHSFxHJIEr6IiIZRElfRCSDKOmLiGQQJX0RkQyipC8ikkGU9EVEMoiSvohIBlHSFxHJIEm3iIqZHQTej9Ph2gAfx+lY9UUxxkeyx5js8YFijBe/YrzIOdc2WqGkS/rxZGYbvKwk4yfFGB/JHmOyxweKMV6SPUY174iIZBAlfRGRDJLuSX+u3wF4oBjjI9ljTPb4QDHGS1LHmNZt+iIiUlG61/RFRCRMxiR9M/uFmTkza+N3LJWZ2Uwz+4eZbTGzF82sld8xAZjZcDPbaWa7zazQ73gqM7OOZrbKzN4xs+1m9q9+x1QdM8sys81m9he/Y4nEzFqZ2fOB7+E7ZnaZ3zGFM7OfBv6Pt5nZM2bW1O+YAMzsSTP7yMy2hW07z8xeNrNdgX/P9TPGyjIi6ZtZR2AYsNfvWKrxMtDHOdcXeBe4z+d4MLMsYDZwLdALGGtmvfyNqoozwM+dcz2BgcBPkjDGoH8F3vE7iBr8FljunOsBXEoSxWpm7YF7gHznXB8gCxjjb1Qh84HhlbYVAq8657oCrwaeJ42MSPrAY8AvgaTswHDO/dU5dybwdA3Qwc94AgYAu51z7znnvgCeBUb4HFMFzrkPnXObAo8/ozxRtfc3qqrMrANwHTDP71giMbMWwGDgvwGcc1845474G1UVDYFmZtYQOAv4wOd4AHDO/Q04XGnzCGBB4PECYGRCg4oi7ZO+md0A7HfOve13LB79EHjJ7yAoT577wp6XkoQJNcjMsoE8YK2/kUT0OOWVji/9DqQaFwMHgf8daIKaZ2Zn+x1UkHNuP/Abyq/UPwQ+dc791d+oanSBc+5DKK+YAOf7HE8FaZH0zeyVQFtf5Z8RwAPAg0keY7DMA5Q3WSz0L9IQi7AtKa+UzKw58GfgXufcUb/jCWdm1wMfOec2+h1LDRoC/YA/OOfygM9JoiaJQJv4CKAz0A4428zG+xtV6mrodwDx4JwbGmm7meVQ/kV528ygvNlkk5kNcM79vwSGWG2MQWZ2K3A9cLVLjnG0pUDHsOcdSJJL6nBm1ojyhL/QOfeC3/FEMAi4wcy+DTQFWpjZ0865ZEpapUCpcy54lfQ8SZT0gaHAP51zBwHM7AXgcuBpX6Oq3gEzu9A596GZXQh85HdA4dKipl8d59xW59z5zrls51w25V/ufolO+NGY2XDgV8ANzrnjfscTsB7oamadzawx5R1nS3yOqQIrP5P/N/COc+5Rv+OJxDl3n3OuQ+D7NwZYmWQJn8Dfwz4z6x7YdDWww8eQKtsLDDSzswL/51eTRB3NESwBbg08vhX4vz7GUkVa1PTTwCygCfBy4IpkjXNukp8BOefOmNldwArKR0s86Zzb7mdMEQwCvg9sNbPiwLb7nXPLfIwpVd0NLAyc4N8DfuBzPCHOubVm9jywifLmz80kyV2vZvYMcCXQxsxKganADOBPZnY75SesUf5FWJXuyBURySBp3bwjIiIVKemLiGQQJX0RkQyipC8ikkGU9EVEMoiSvohIBlHSFxHJIEr6IiIZ5P8DV6+ZLRYnN1MAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "### M-H采样\n",
    "# 导入相关库 \n",
    "import random\n",
    "from scipy.stats import norm\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# 定义平稳分布为正态分布\n",
    "def smooth_dist(theta):\n",
    "    '''\n",
    "    输入：\n",
    "    thetas：数组\n",
    "    输出：\n",
    "    y：正态分布概率密度函数\n",
    "    '''\n",
    "    y = norm.pdf(theta, loc=3, scale=2)\n",
    "    return y\n",
    "\n",
    "# 定义M-H采样函数\n",
    "def MH_sample(T, sigma):\n",
    "    '''\n",
    "    输入：\n",
    "    T：采样序列长度\n",
    "    sigma：生成随机序列的尺度参数\n",
    "    输出：\n",
    "    pi：经M-H采样后的序列\n",
    "    '''\n",
    "    # 初始分布\n",
    "    pi = [0 for i in range(T)]\n",
    "    t = 0\n",
    "    while t < T-1:\n",
    "        t = t + 1\n",
    "        # 状态转移进行随机抽样\n",
    "        pi_star = norm.rvs(loc=pi[t-1], scale=sigma, size=1, random_state=None)  \n",
    "        alpha = min(1, (smooth_dist(pi_star[0]) / smooth_dist(pi[t-1])))   \n",
    "        # 从均匀分布中随机抽取一个数u\n",
    "        u = random.uniform(0, 1)\n",
    "        # 拒绝-接受采样\n",
    "        if u < alpha:\n",
    "            pi[t] = pi_star[0]\n",
    "        else:\n",
    "            pi[t] = pi[t-1]\n",
    "    return pi\n",
    "    \n",
    "# 执行MH采样\n",
    "pi = MH_sample(10000, 1)\n",
    "\n",
    "### 绘制采样分布\n",
    "# 绘制目标分布散点图\n",
    "plt.scatter(pi, norm.pdf(pi, loc=3, scale=2), label='Target Distribution')\n",
    "# 绘制采样分布直方图\n",
    "plt.hist(pi, \n",
    "         100, \n",
    "         normed=1, \n",
    "         facecolor='red', \n",
    "         alpha=0.6, \n",
    "         label='Samples Distribution')\n",
    "plt.legend()\n",
    "plt.show();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAEICAYAAACktLTqAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAGOZJREFUeJzt3X2UXHV9x/H3p5uELRCBs4kPZIMbJT4E5AjdBEq1RR7kQU3EagmERkVOBEpbtbbYAmkAT8VHhEpVipyAIlFatakGoRaop9pglgeBBKghBlkCMYQSEiBIyLd/3Lt4mczs3t2dp/3t53XOntw79zf3fuc3k8/c+d07dxQRmJlZWn6n1QWYmVn9OdzNzBLkcDczS5DD3cwsQQ53M7MEOdzNzBLkcLdhkbRE0jfy6f0kbZPU0YDtLJX0yUGWb5P0mnpvdyyRdKakjXlfdLW6ngGSbpV0eqvrGO8c7mOEpLdI+qmkLZKekPQTSbNbWVNE/Coi9oyIF1qw7T0jYt1gbSQdIam/WTU1k6SJwBeAt+d9sbnVNVl7mdDqAmxokl4GfB84E/g2MAl4K/BcK+sa7yR1tOKNLfcKoBNY3aLtW5vznvvY8DqAiLguIl6IiGcj4qaIuBtA0msl3Sxps6THJV0rae+BO0taL+mvJd0t6WlJX5P0Ckk3SNoq6UeS9snb9kgKSYskbZD0qKS/qlZUoe2EfP5WSRflnyq2SrpJ0pRC+4WSHsrrPD+v6+hBHvc+kn6Qr+s2Sa8trCsk7Z9PnyBpTd7uEUkfl7QHcAOwbz5ssU3SvpJ2k/TF/LFtyKd3K6z3b/LHvEHS6RXbWSrpy5JWSHoaeJukd0i6U9JTkh6WtKRK/3wwX/Z/ks6QNDt/Lp6U9KVaD75WrZJeBzyQN3tS0s1V7nuSpHX5jgGSjpf0mKSpNbZ1fb58i6QfSzqgsGyppMsHeS6OkXR/ft8vAaqxjVdKekaFISRJvydpU/5JxOopIvzX5n/Ay4DNwNXA8cA+Fcv3B44BdgOmAj8GvlhYvh5YSba3Nw34NXAHcHB+n5uBv8/b9gABXAfsAbwJ2AQcnS9fAnyjou2EfP5W4EGyN6PfzecvzpfNArYBbyH75PE54PmB9VZ5zEuBJ4A5ZJ8wrwWWFZYHsH8+/Sjw1nx6H+CQfPoIoL9ivRfmffHyvK9+ClyULzsOeAw4ANgd+HrFdpYCW4A/INsx6sy38aZ8/iBgI/Duiv75St727cB24Hv59geeiz+q0QeD1fqSvq9x/2vzmruADcA7B2l7GjA5fz18EbirzHMBTAGeAt4LTAQ+CuwATq+xnRXAmYX5S4B/bPX/sRT/Wl6A/0o+UfDG/D9Zf/6fZznwihpt3w3cWZhfDywozP8r8OXC/J8D38unB0LjDYXlnwG+lk8vYfBwP69wv7OAH+bTi4HrCst2B37D4OF+ZWH+BOD+wnwxdH8FfBh4WcU6jmDXcH8QOKEwfyywPp++CvhUYdn+7Bru1wzxPH0RuKSif6YVlm8GTqp4Lj5SY12D1fqSvq9x/73zvrkH+OowXmt75+vea6jnAlgIrCwsE9lrtFa4nwT8JJ/uIHszndPq/18p/nlYZoyIiPsi4gMR0Q0cCOxLFiRIermkZfmQxFPAN8j2qIo2FqafrTK/Z0X7hwvTD+XbK+OxwvQzhfXuW1xnRDxDFnQjWVelPyYLnIck/Zek3x9knfuSPZ4Bxcf2khorpqveJulQSbfkQwtbgDMYfd+XqXVIEfEkcD3Z6+XztdpJ6pB0saQH89fP+nxR8XGUfV6D6v024N+AWcrOdDoG2BIRPyv3iGw4HO5jUETcT7Y3dWB+06fI9rQOioiXAadSY9xzGKYXpvcj+1g/Go8C3QMzkn6XbLhg1CJiVUTMIxu++B7ZQWfI+qTSBuDVhfniY3tJjby0D17cXMX8N8k+RU2PiL3IhmBG2/cDBqt1SJLeTDbcch1w2SBNTwHmAUcDe5F9KoByj+NRCv0kSVTvNwAiYjvZ87MA+FOyoS9rAIf7GCDpDZL+SlJ3Pj8dOJlsPBaysdJtZAfXpgF/XYfNni9p9/zA2geBb41yff8CvEvS4ZImARdQhxCUNEnSAkl7RcTzZOO/A2ewbAS6JO1VuMt1wHmSpuYHexeTfdKBLHQ+KOmNknbPlw1lMvBERGyXNIcsKOtlsFoHJakzb/t3ZM/fNEln1Wg+mezMq81kw2X/MIwafwAcIOk9yg6s/wXwyiHucw3wAWAuJR+PDZ/DfWzYChwK3JafpbESuBcYOIvlAuAQsoN9PwC+U4dt/hewFvhP4HMRcdNoVhYRq8nG9peR7e1tJTuYWI/TOf8UWJ8PKZxB9sll4BPOdcC6/MyUfYFPAn3A3WRj0XfktxERN5Dt4d5C9tj/J1//YDWeBVwoaStZ+H57kLbDVbPWEj5FdrzhyxHxHFmffFLSzCptryEb8nkEWMNvdxqGFBGPA+8DLiZ7c5gJ/GSI+/wE2AncERHry27Lhkf5gQ0zIDt9D/glMDEidjRwO3sCTwIzI+KXjdrOaEh6I9mb6G6N7IvxKD9985sRcWWra0mV99ytaSS9Kx/q2YPsVMh7+O3Bu7Yg6cR8qGcf4NPAvzvY60vZN6sPYfRDfTYIh7s10zyyA4IbyD6+z4/2++j4YbLz+h8kG7s/s7XlpEXS1cCPyE7/3NrqelLmYRkzswR5z93MLEEtu3DYlClToqenp1WbNzMbk26//fbHI6LqNYKKWhbuPT099PX1tWrzZmZjkqSHhm7lYRkzsyQ53M3MEuRwNzNLkH+JyczGneeff57+/n62b9/e6lJq6uzspLu7m4kTR/Y7Jg53Mxt3+vv7mTx5Mj09PWQXsmwvEcHmzZvp7+9nxowZI1qHh2XMbNzZvn07XV1dbRnsAJLo6uoa1ScLh7uZjUvtGuwDRlufw93MLEEeczczW7Kkvdc3At5zN7PGW7KkLQJvPHG4m5k12apVqzjooIPYvn07Tz/9NAcccAD33ntvXbfhYRkza57i3vs43pOfPXs2c+fO5bzzzuPZZ5/l1FNP5cADDxz6jsPgcDcza4HFixcze/ZsOjs7ueyyy+q+fg/LmJm1wBNPPMG2bdvYunVrQ74p63A3M2uBRYsWcdFFF7FgwQLOOeecuq/fwzJmZk0e/7/mmmuYMGECp5xyCi+88AKHH344N998M0ceeWTdtuFwN7PWGMcHVxcuXMjChQsB6Ojo4Lbbbqv7NhzuZtYY4yyw243H3M3MEuRwNzNLkMPdzCxBDnczswQ53M3MEuSzZcxs3LtkZX3X99HD6ru+kfCeu5lZghzuZmZNdv7553PppZe+OH/uuefW/eJhDnczsyb70Ic+xNVXXw3Azp07WbZsGQsWLKjrNjzmbmbWZD09PXR1dXHnnXeyceNGDj74YLq6uuq6DYe7mdWPLzlQ2umnn87SpUt57LHHOO200+q+fg/LmJm1wIknnsgPf/hDVq1axbHHHlv39XvP3czGvVacujhp0iTe9ra3sffee9PR0VH39Tvczaz1BoZzxtGwzs6dO1m5ciXXX399Q9bvYRkzsyZbs2YN+++/P0cddRQzZ85syDZK7blLOg64FOgAroyIi2u0ey9wPTA7IvrqVqWZWUJmzZrFunXrGrqNIffcJXUAlwPHA7OAkyXNqtJuMvAXQP1/UsTMrM4iotUlDGq09ZUZlpkDrI2IdRHxG2AZMK9Ku4uAzwD1/xlvM7M66uzsZPPmzW0b8BHB5s2b6ezsHPE6ygzLTAMeLsz3A4cWG0g6GJgeEd+X9PFaK5K0CFgEsN9++w2/WjOzOuju7qa/v59Nmza1upSaOjs76e7uHvH9y4S7qtz24tudpN8BLgE+MNSKIuIK4AqA3t7e9nzLNLPkTZw4kRkzZrS6jIYqMyzTD0wvzHcDGwrzk4EDgVslrQcOA5ZL6q1XkWZmNjxlwn0VMFPSDEmTgPnA8oGFEbElIqZERE9E9AArgbk+W8bMrHWGDPeI2AGcDdwI3Ad8OyJWS7pQ0txGF2hmZsNX6jz3iFgBrKi4bXGNtkeMviwzMxsNf0PVzCxBDnczswQ53M3MEuRwNzNLkMPdzCxBvp67mY3OOLoG+1jiPXczswQ53M3MEuRwNzNLkMPdzCxBPqBqZu2jeHDWB2pHxXvuZmYJcribmSXI4W5mliCHu5lZghzuZmYJcribmSXI4W5mliCHu5lZghzuZmYJcribmSXI4W5mliCHu5lZghzuZmYJcribmSXI4W5mliCHu5lZghzuZmYJ8i8xmdnI+JeS2prD3czak39yb1Q8LGNmliCHu5lZghzuZmYJcribmSXI4W5mliCHu5lZgkqFu6TjJD0gaa2kT1RZfoakeyTdJem/Jc2qf6lmZlbWkOEuqQO4HDgemAWcXCW8vxkRb4qINwOfAb5Q90rNzKy0Mnvuc4C1EbEuIn4DLAPmFRtExFOF2T2AqF+JZmY2XGW+oToNeLgw3w8cWtlI0p8BHwMmAUdWW5GkRcAigP3222+4tZqZWUll9txV5bZd9swj4vKIeC1wDnBetRVFxBUR0RsRvVOnTh1epWZmVlqZcO8Hphfmu4ENg7RfBrx7NEWZmdnolAn3VcBMSTMkTQLmA8uLDSTNLMy+A/hF/Uo0M7PhGnLMPSJ2SDobuBHoAK6KiNWSLgT6ImI5cLako4Hngf8D3t/Ios3MbHClLvkbESuAFRW3LS5M/2Wd6zIzs1HwN1TNzBLkcDczS5DD3cwsQQ53M7MEOdzNzBLkcDczS5DD3cwsQQ53M7MElfoSk5kZAEuWtLoCK8nhbmbtr/im4jeYUjwsY2aWIIe7mVmCHO5mZglyuJuZJcjhbmaWIIe7mVmCHO5mZglyuJuZJcjhbmaWIIe7mVmCHO5mZglyuJuZJcjhbmaWIIe7mVmCHO5mZglyuJuZJcjhbmaWIIe7mVmCHO5mZglyuJuZJcg/kG1mg2u3H6QeqKfd6moz3nM3M0uQw93MLEEOdzOzBDnczcwS5HA3M0tQqbNlJB0HXAp0AFdGxMUVyz8GnA7sADYBp0XEQ3Wu1ayUS1buettHD2t+HWatNGS4S+oALgeOAfqBVZKWR8SaQrM7gd6IeEbSmcBngJMaUbDZSFQLfHDoW7rK7LnPAdZGxDoAScuAecCL4R4RtxTarwROrWeRZtXUCmwzKzfmPg14uDDfn99Wy4eAG6otkLRIUp+kvk2bNpWv0szMhqVMuKvKbVG1oXQq0At8ttryiLgiInojonfq1KnlqzQzs2EpMyzTD0wvzHcDGyobSToaOBf4o4h4rj7lmZnZSJQJ91XATEkzgEeA+cApxQaSDga+ChwXEb+ue5VmDeIDrZaqIcM9InZIOhu4kexUyKsiYrWkC4G+iFhONgyzJ3C9JIBfRcTcBtZt44gPnJoNX6nz3CNiBbCi4rbFhemj61yXmZmNgr+hamaWIIe7mVmCHO5mZgnyLzGZVeGzaGys8567mVmCvOdubcOnPJrVj8PdzHblH58e8zwsY2aWIIe7mVmCPCxjZmNTcejIw0i78J67mVmCHO5mZgnysIzZMPjLTTZWeM/dzCxB3nO3pvOXlcwaz3vuZmYJcribmSXI4W5mliCHu5lZghzuZmYJcribmSXIp0Ka1YG/3GTtxnvuZmYJ8p67NYy/rGTWOt5zNzNLkMPdzCxBHpYxs98aqz964R/u2IX33M3MEuRwNzNLkMPdzCxBDnczswQ53M3MEuSzZWzU/GWl2nxZAmsV77mbmSXI4W5mliCHu5lZgkqFu6TjJD0gaa2kT1RZ/oeS7pC0Q9J761+mmZkNx5DhLqkDuBw4HpgFnCxpVkWzXwEfAL5Z7wLNzGz4ypwtMwdYGxHrACQtA+YBawYaRMT6fNnOBtRoZmbDVGZYZhrwcGG+P79t2CQtktQnqW/Tpk0jWYWZmZVQJtxV5bYYycYi4oqI6I2I3qlTp45kFWZmVkKZYZl+YHphvhvY0JhyzKzpfIncJJUJ91XATEkzgEeA+cApDa3KLHH+5moDDbxZjfM3rSHDPSJ2SDobuBHoAK6KiNWSLgT6ImK5pNnAd4F9gHdJuiAiDmho5dZ0vsyA2dhR6toyEbECWFFx2+LC9Cqy4RozM2sD/oaqmVmCHO5mZglyuJuZJcjhbmaWIIe7mVmCHO5mZgnyz+yZtZFq3yXwF5tsJBzuVpW/sGQ2tnlYxswsQd5zNxuPxvl1V8YDh7uZpan4BjYO38w8LGNmliCHu5lZghzuZmYJ8pi7WZvzD3vYSDjcxzmfz26WJg/LmJklyOFuZpYgD8uYjSfj8Hzv8crhbmbpG4dfaPKwjJlZgrznPk74rJj0+BRJG4z33M3MEuRwNzNLkMPdzCxBDnczswQpIlqy4d7e3ujr62vJtlPmA6e2i1tvBeCj229taRltaQyeFinp9ojoHaqd99zNzBLkcDczS5DD3cwsQf4S0xjlsXUbrks6j6h6+7geix8Ycx+DY+9DcbibpSg/iGrjl4dlzMwS5D13s3HOwzVpcri3OY+tW2keihm5BC8J7HBvEw5xazfeox/bSoW7pOOAS4EO4MqIuLhi+W7ANcDvAZuBkyJifX1LNbOqvMdeX4nsxQ95+QFJHcD/AscA/cAq4OSIWFNocxZwUEScIWk+cGJEnDTYesfr5Qe8h251MYYCPZk9/TYJ+rKXHyiz5z4HWBsR6/IVLwPmAWsKbeYBS/LpfwG+JEnRqgvXNIiD2ZpuDIV4LbWGd6pp6zeCauHeJoFfTZlwnwY8XJjvBw6t1SYidkjaAnQBjxcbSVoELMpnt0l6IJ+eUtm2TbnO+nKd9TMWaoQh6vxYEwsZQrn+vOCCxleyq1eXaVQm3FXltso98jJtiIgrgCt22YDUV+ZjRqu5zvpynfUzFmoE19lMZb7E1A9ML8x3AxtqtZE0AdgLeKIeBZqZ2fCVCfdVwExJMyRNAuYDyyvaLAfen0+/F7g5tfF2M7OxZMhhmXwM/WzgRrJTIa+KiNWSLgT6ImI58DXg65LWku2xzx9mHbsM1bQp11lfrrN+xkKN4DqbpmW/xGRmZo3jC4eZmSXI4W5mlqCmhbuk90laLWmnpN6KZX8raa2kByQdW+P+MyTdJukXkr6VH9xtdM3fknRX/rde0l012q2XdE/erulfu5W0RNIjhVpPqNHuuLyP10r6RAvq/Kyk+yXdLem7kvau0a7p/TlU30jaLX89rM1fhz3NqKuihumSbpF0X/5/6S+rtDlC0pbCa2Fxs+vM6xj0OVTmsrw/75Z0SAtqfH2hn+6S9JSkj1S0aYv+HJGIaMof8Ebg9cCtQG/h9lnAz4HdgBnAg0BHlft/G5ifT38FOLNZtefb/DywuMay9cCUZtZTsf0lwMeHaNOR9+1rgEl5n89qcp1vBybk058GPt0O/Vmmb4CzgK/k0/OBb7XgeX4VcEg+PZnssiCVdR4BfL/ZtQ33OQROAG4g+47MYcBtLa63A3gMeHU79udI/pq25x4R90XEA1UWzQOWRcRzEfFLYC3ZJQ9eJEnAkWSXNgC4Gnh3I+utsv0/Aa5r1jYb4MXLSETEb4CBy0g0TUTcFBE78tmVZN+ZaAdl+mYe2esOstfhUfnromki4tGIuCOf3grcR/bt8LFoHnBNZFYCe0t6VQvrOQp4MCIeamENddUOY+7VLm9Q+YLtAp4sBEO1No30VmBjRPyixvIAbpJ0e36JhVY4O/94e5WkfaosL9PPzXQa2Z5bNc3uzzJ985JLbAADl9hoiXxY6GDgtiqLf1/SzyXdIOmAphb2W0M9h+32epxP7Z23dujPYavr9dwl/Qh4ZZVF50bEv9W6W5XbRnR5g5EoWfPJDL7X/gcRsUHSy4H/kHR/RPy4HvWVqRP4MnARWZ9cRDaEdFrlKqrct+7nwZbpT0nnAjuAa2uspuH9WaGlr8HhkrQn8K/ARyLiqYrFd5ANLWzLj718D5jZ7BoZ+jlsp/6cBMwF/rbK4nbpz2Gra7hHxNEjuFuZyxs8TvaxbUK+11StzYgMVbOyyym8h+xa9bXWsSH/99eSvkv2Mb+uYVS2byX9M/D9KovK9POolejP9wPvBI6KfFCzyjoa3p8VhnOJjX618BIbkiaSBfu1EfGdyuXFsI+IFZL+SdKUiGjqRcVKPIdNeT2WdDxwR0RsrFzQLv05Eu0wLLMcmJ+fjTCD7F3xZ8UGeQjcQnZpA8gudVDrk0C9HQ3cHxH91RZK2kPS5IFpsoOG9zaptoEaimOVJ9bYfpnLSDSUsh99OQeYGxHP1GjTiv4cE5fYyMf4vwbcFxFfqNHmlQPHAiTNIfs/vrl5VZZ+DpcDC/OzZg4DtkTEo82ss6DmJ/N26M8Ra9aRW7LQ6QeeAzYCNxaWnUt2tsIDwPGF21cA++bTryEL/bXA9cBuTap7KXBGxW37AisKdf08/1tNNvzQ1KPiwNeBe4C7yf7TvKqyznz+BLIzLB5sUZ1rycZZ78r/vlJZZ6v6s1rfABeSvREBdOavu7X56/A1Lei/t5ANXdxd6MMTgDMGXqPA2Xm//ZzsoPXhLaiz6nNYUaeAy/P+vofCGXRNrnV3srDeq3BbW/XnSP98+QEzswS1w7CMmZnVmcPdzCxBDnczswQ53M3MEuRwNzNLkMPdzCxBDnczswT9Pw9SK358vN6hAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "### Gibbs采样\n",
    "# 导入math库\n",
    "import math\n",
    "# 导入多元正态分布函数\n",
    "from scipy.stats import multivariate_normal\n",
    "\n",
    "# 指定二元正态分布均值和协方差矩阵\n",
    "target_distribution = multivariate_normal(mean=[5,-1], cov=[[1,0.5],[0.5,2]])\n",
    "\n",
    "# 定义给定x的条件下y的条件状态转移分布\n",
    "def p_yx(x, mu1, mu2, sigma1, sigma2, rho):\n",
    "    '''\n",
    "    输入：\n",
    "    x：公式(25.30)中的x2\n",
    "    mu1：二维正态分布中的均值1\n",
    "    mu2：二维正态分布中的均值2\n",
    "    sigma1：二维正态分布中的标准差1\n",
    "    sigma2：二维正态分布中的标准差2\n",
    "    rho：公式(25.30)中的rho\n",
    "    输出：\n",
    "    给定x的条件下y的条件状态转移分布\n",
    "    '''\n",
    "    return (random.normalvariate(mu2 + rho * sigma2 / sigma1 * (x - mu1), math.sqrt(1 - rho ** 2) * sigma2))\n",
    "\n",
    "# 定义给定y的条件下x的条件状态转移分布\n",
    "def p_xy(y, mu1, mu2, sigma1, sigma2, rho):\n",
    "    '''\n",
    "    输入：\n",
    "    y：公式(25.31)中的x1\n",
    "    mu1：二维正态分布中的均值1\n",
    "    mu2：二维正态分布中的均值2\n",
    "    sigma1：二维正态分布中的标准差1\n",
    "    sigma2：二维正态分布中的标准差2\n",
    "    rho：公式(25.31)中的rho\n",
    "    输出：\n",
    "    给定y的条件下x的条件状态转移分布\n",
    "    '''\n",
    "    return (random.normalvariate(mu1 + rho * sigma1 / sigma2 * (y - mu2), math.sqrt(1 - rho ** 2) * sigma1))\n",
    "\n",
    "def Gibbs_sample(N, K):\n",
    "    '''\n",
    "    输入：\n",
    "    N：采样序列长度\n",
    "    K：状态转移次数\n",
    "    输出：\n",
    "    x_res：Gibbs采样x\n",
    "    y_res：Gibbs采样y\n",
    "    z_res：Gibbs采样z\n",
    "    '''\n",
    "    x_res = []\n",
    "    y_res = []\n",
    "    z_res = []\n",
    "    # 遍历迭代\n",
    "    for i in range(N):\n",
    "        for j in range(K):\n",
    "            # y给定得到x的采样\n",
    "            x = p_xy(-1, 5, -1, 1, 2, 0.5)  \n",
    "            # x给定得到y的采样\n",
    "            y = p_yx(x, 5, -1, 1, 2, 0.5)   \n",
    "            z = target_distribution.pdf([x,y])\n",
    "            x_res.append(x)\n",
    "            y_res.append(y)\n",
    "            z_res.append(z)\n",
    "    return x_res, y_res, z_res\n",
    "\n",
    "# 二维正态的Gibbs抽样           \n",
    "x_res, y_res, z_res = Gibbs_sample(10000, 50)\n",
    "\n",
    "# 绘图\n",
    "num_bins = 50\n",
    "plt.hist(x_res, num_bins, normed=1, facecolor='red', alpha=0.5, label='x')\n",
    "plt.hist(y_res, num_bins, normed=1, facecolor='dodgerblue', alpha=0.5, label='y')\n",
    "plt.title('Sampling histogram of x and y')\n",
    "plt.legend()\n",
    "plt.show();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAb4AAAEuCAYAAADx63eqAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzsvXmQHNd95/mtu6rvRl9odDWABhpXd6MPEqBJOzyWRWnggRSIYAwdpiIo2qJoKxTUBmk7wqIvUmJo12tbM94Ic//w7GpM2msuaGp3BiMflMi1LFq8QFKNq+8T6Puo+8iqymv/aL5kVlZmVmYd2dVV7xPRga6urHxVhcz3fb/f+x02URRBoVAoFEqtYN/vN0ChUCgUipVQ4aNQKBRKTUGFj0KhUCg1BRU+CoVCodQUVPgoFAqFUlNQ4aNQKBRKTeHM8zzNdaBQKBTKQcSm9QS1+CgUCoVSU1Dho1AoFEpNQYWPQqFQKDUFFT4KhUKh1BRU+CgUCoVSU1Dho1AoFEpNQYWPQqFQKDUFFT4KhUKh1BRU+CgUCoVSU1Dho1AoFEpNQYWPQqFQKDUFFT4KhUKh1BRU+CgUCoVSU1Dho1AoFEpNQYWPQqFQKDUFFT4KhUKh1BRU+CgUCoVSU1Dho1AoFEpNQYWPQqFQKDUFFT7KviCK4n6/BQqFUqM49/sNUGqPTCaDTCYDp9MJm82m+gMg518KhUIpBVT4KJYhiiJSqRTS6TQcDof0N/IDAIuLi2hvb0dTUxOAbPGz2/ccFORfNZGkUCiUfFDho1iCIAhIJBLgOA4OhwN2uz1LrMjvHMeB53lJ3IggiqIIjuM0z69lOcqFkYojhUIBqPBRLCCTySCZTAIAHA6HrgDZbLas/T81cVQiF0e59ah2bvKv3HqkrlUKpbagwkcpG6IogmEYpNNp2O12SWz0sNvtEATB1DhKoVITLrkYEusxGo1ifX0dZ8+ezXm9/Ie6VimU6oIKH6Us8DyPZDIpuTaNCkUhwmcENXG02WwQBCFLkOXWo977UHOnUuuRQjkYUOGjlBRRFMGyLObm5tDY2IjW1lZTr7fb7ZalOijdquRvar/LMepaJecgVqMoipI4UuuRQtk/qPBRSobctcmyLHieN30OYoVZgZrwGX2d3mOCMjDn1q1bOHnyJOrr63NeTwNzKBTroMJHKQk8zyORSEAQBDgcDjgcDl0B0xKccrk61ShU+MycX/67cs8QoIE5FMp+QIWPUhSiKCKTyYBhGACQ8vP0BEwURayvr4NlWXg8Hni9XumHPG8F5RY+o+9B7zGgHpijdz7yw/M83G531t+1xqBQagkqfJSCIa7NTCajmpenJnwsy2Jqago+nw8dHR1IpVKIRCLY2tpCKpVCIpGAKIrY2NjIEUWv1wuPxwOPx2MoQjQflSB8RijEtSoIAt577z089NBDmueT/59R65FSS1DhoxSE3LWpFD1A3eKLRCKYmZlBX18fOjs74XK5cs67sbGBVCqF3t5eqcoLEcRAIIBUKoVMJgNRFGGz2XLEUf6YWJ9aWC18++FaVS4Q5OMLgmA4MAegFXMo1QMVPoop5K5Nm82mKS52u11yyYmiiJWVFezs7OD8+fPw+Xya5yeWotPpRENDAxoaGjSPFQQhSxyJ9UgeE+F1u905okjcqlbtJ1YKhQbm6J2P/HAcB7fbTQNzKBUPFT6KYQRBQDKZxMTEBM6dO6c7oRGLj7g2vV4vxsbG8roozQS32O121NXVoa6uTvMYItREGNPpNEKhEFKpFJLJJOLxON555x24XC5VcSSu1VqbvAtJ63j//fd1Xas0MIdSKVDhoxiC4zgkk0kIgoB4PJ53krLb7WAYBuPj4+jr60NHR0fW83rRi6W0wog71OPxoLm5Oeu5TCaDmzdv4sKFC+A4LkscY7EYdnZ2kE6nkU6npc8kF0alSJZi3/EgoWY96rlWzdZbpa5VSrmgwkfRxahrU/maQCCAQCCA++67T9e1qcTqBHbyr8vlgsvlQmNjo+bxPM9nuVVTqRTC4bD0N/K+tcSx1tyqQOGBOfnOR1zpTqczqzIQFUiKEajwUTQhrk2WZQ2XHWNZFtPT0xBFER0dHaZED7A+j8/MWA6HI69rVRCELNdqKpWSgnKCwSBisRgcDgdcLpdmYA7ZJ6slzLhWgb3veXx8HIODg/B4PKrno4E5FC2o8FFUYVkWyWQSoigaFr1oNIrp6WkcP34cbrcb29vbuseTyEw5Vlp85cBut2cFz8i5ffs2jh07hsbGRrAsmxWYE4lEJKFkWRbAntBqiWOpUjoOEsprhQRBablXaSsrihZU+ChZiKKIdDoNhmFgt9sNuzZXV1exvb0tRW1Go9GCLDcrS5btl8jabDa43W4puVwLsu8oT+kIBoPS35QpHUqRPMgLCCMoC4wTiq23Kk87UVqHoijC5XJR1+oBhwofRaIY16bH48mK2izUZWmlqxOwrkoMYH5yNJrSId93JIE5JGr1nXfekYRWK2pVPpEXyn6IrJrHwCj59h6VOZeiKGJrawvRaBT9/f2q56OBOQcHKnwUAHsCtrCwgEOHDsHr9erepMQqi8fjkmuzs7Mz6xgjAqY2hpVJ5ftRuaXU49ntdvh8PtW91HfffRc///M/n5XSQUSSpHSQguKiKMLpdGoWA8iX0lGMCBUK6XZRLpSfh3xHeq5Vo4E5tGLO/kKFr8aRuzYjkQiampry3ng2my1vQroR4dPa47MyuKUWkKd0aEH2w+TWYzwex+7ublZKh81mUxVHtSo81QbP86qu/0ICc8y0siJjUuuxdFDhq2GIa5M0i83XUQHYswwZhkEymcR9992nueIuxtVZ7XtTlYg8pUPPtapM6Uin04hEImAYRioGIBdaNZF0Op0HcuImnUcKodCKOTdu3MCZM2dyFpc0MKc4qPDVKCQ3T94cNZ9YkahNr9eLvr4+XTdToUEqVga3UMyjldLBsizGx8fxwAMPSF4EuWuVpHSkUimpT+NBS+ngeR5OZ3mnTKX1yPM8XC6X4VZWWoE5NK0jGyp8NYYoitKEpExIt9vtqs1jSdTm1tYWhoaGsLi4mFecCrXcrA5usZL92FO0CrnbWu4O1TueZVnNOqukELnD4dCslGP1d1mMxVcoamJbSGCOVloHCXyqtdQYKnw1hCAISCQSkmtTecOouTo5jsPU1BQ8Ho/k2jQiTgdhj49SOrRSC7SQp3ToVctR7jsyDCMF5iQSCcm1qhWQQ1yrpYDnecsFQhCEgiwyo67VWr3XqPDVCJlMBslkEgA0UxWUoiNPSJdHbRqN2My3It/vqE5K6ShXVKfT6YTT6UR9fX3OcyRqVZ7SQf4lKR3pdBocx0l7mFpRq0ZSOrSCW8pNLboiyw0VvipHFEXs7u6C53n4fD7dFSsRNFEUsba2hs3NTQwNDeXs55TTKqM3+cHE6nQG+eJIL6VDfnwmk8myHsPhsG5Kh9K1Sjwl1UYt3nNU+KoYnueRTCaxubkJj8ejumqW43A4kE6ncefOHbjdboyNjane6EaFrxZvKD2q2Zq1WvgKca2SSNOmpibN4+RdOpQNkMPhMMLhsCSOWoE51SiO1QYVviqEBA4Q16bL5TIkVOl0Gqurqzh9+nROQrqcUlh8oigiHA5Lk0Wtba6Xm/1IzLfy/9Cs8BlFr1rO+Pg4zpw5A6/Xm5XOIa+1mk6nwfN8VrUcNQvSSEpHuf8Pq3URZgQqfFWGKIpgGAbpdFray5N3Q9d6zdraGjY2NtDZ2akrekDxwseyLCYnJ6WgA9LSh6zKU6kUFhcXsyYNKo6VzX64Ovcj0MThcBhugEwS/4lIBoNB6TEpRC7fd1Tr72hF+kQtemao8FURPM8jkUhINyi5oIkLUw2O4zA9PQ2Xy4VTp04hHA7nHacY4SMBMydOnMDhw4ezJi9BEJBKpfDRRx/B5/PldC1QFmVW/ng8Hupm+oT9cD1W83iAuahOeaSpsgEygaQZyF2rpAGy/G/vv/++rmuVLgjNQ4WvCpA3iwWQM/lrCVUsFsPU1BSOHTuGrq4uhMNhQ4JmRvjIBCwPmNErc1ZXVwen04nu7m7V8ymLMquJIwBNcdzPPZhqXllX+h5fKSh1VKe8Wo5aSkc8Hsf8/DzOnz+f5VqVp3Sk02lpEaAljKVM6agW6LdxwCGuzUwmk1XwVo7D4chKTNeK2jQqaGaCW0jh3unpadjtds2AGaMYieBTawZLwttTqZT03j0eDxiGwcLCQs6EUa6JopK7QRRDLQif1e5VkrzucDhQX1+vG5xGFoTyRWG+lA6Px4Ouri7dQgPVChW+A4zctaklekC2UMldm0oRUgqkFvn2DAk2mw3xeBwzMzPo6enBkSNHDH6y4tBrBksgezDXr19HfX29rjhqWY+VvIrej+CWahc+qzFjYRpN6SANkJUeklqjcu9ciibEtbmwsAC3242uri7d44mgKV2bSoxackYFklR9GRgY0K3OoUa5J1KyB+N0OtHV1aU6llqAglrNSWX0XiWIY7ULUS0IH8dxJb1+5NVympqa9qUEW6VAhe+AIQiC5NokRWzzIbe81BLSCaVydQqCgIWFBWQyGdx3332aLhqt1Sap87nfe2JGAxSUblWtgswkmi+RSOSI435/1mKpdqHdD6yoFHPQr7tCocJ3gOA4DslkUlqpORyOvC5HjuMwPz8PlmXx4IMP6t5IpRC+dDqNiYkJtLW1obm5uaDJiXRosGJiI/uQhU4A8sRoPXFkWRZTU1Ooq6uDKIpScEIqlZL+D5UltYoVR7rHd7DZrxJptQAVvgOAPGpT3lFBL00B2IvanJ6eRk9Pj6GbyMwen5rwBYNBzM/P49SpU2htbcXExETBPfmsbEZb7n0O4mLyeDxobW1FW1tbzjFq+y9a4qi15yivN0n3+ErLfuyFcRxXEw1+9wMqfBUOaRbLsmxOcWktgRBFEevr69jY2MDg4CC8Xi/W19fzjmVUBJTjiqKIu3fvIhgMYmRkROr0fRCa0VZKGTHl/osaylY+DMMgEolga2tLcqOKogiXywW32w2GYXDv3j1NcSwl1S58+9WSqNwRl9VuNWtBha+CIWXHSF8ytTZCSguNRG06nc6sqM1S1taUC5rchTc6Opp1IxUqfFY2o60U4TOCkVY+JCmaBDLZbLYccQSQU2/S5/MV1QS22oNb9qMlUamDWyifQr/VCoREEzIMA7vdrrnSVAofcW329vbi8OHDZXt/RNDI5NrX14eOjo6c4woVsGpzdVo5FsnVamhogMfjQW9vr+pxcrdqOp3OaQILIKsJrPJHKY7VbvHtx34b3eMrH1T4KgxBELCxsSHtCelNJkT45K7NgYEB1SjKUk5KNpsNyWTSUJSo3kSfL6rTCg6SxWeGfJ9Jr2IIQa+cllwcPR6P1LJHWQigXGJY7RYmoN59vdTQqE7KvkNcm6urq+jt7c3r37fb7VLBZ4fDUXRVFCPwPI+FhQWwLIsLFy6UJEpUCXV1loZiJzW9TgUEIo4rKytSqgZJ5yCBV/KCAmr1VQt5n4IgWOoG3A/rq1r7/1UCVPgqAKVr0+l0Gpr4SXDDmTNnyuraJCSTSUxOTuLw4cNgGCbvTVlMcEu1Cp9VY1k1DhHH+vp62Gw2Vdcqz/OaPe6U4ujxeODz+bL2Hz0eT461tR/lw6irs3qgwrfPkKhNsroj6Qr52gitr69jfX0dPp/PsOgZcddo5dDt7OxgaWkJZ8+eRWNjIzY2NvKOV4zFV42uzmp2KwmCoBl6b6TWpJo4BoNBSRzJHiKxFBOJBOrr67NEspxCuB9RnVYEt1TzNakHFb59hOTmkdWrvI2QVj4dx3GYmZmRXJvj4+OGxnI4HIaEj4gVOU4QBCwuLiKRSGBsbMxUXhG1+PaP/eiPV8x4RsWRBOAsLy+DZVlsbGxoiqPcgiy2hc9+RHXWQpL+fkGFbx8QRVG6WeUJ6QQt4YvH45iamiooarOQqizpdBqTk5NobW3F8PBwQZVDqPDtHwdJ+IzgcDikBrA7Oztoa2tDe3t71jGkp6M8YjUcDkuCma+no5Y47pfbsVzfabVe80ahwmcxgiAgkUhkuTaVqLURyhe1mQ+jVVnIcaFQCHNzc+jv78ehQ4dMjwfk7+KgdfNVq6vTSmq1couR7uhGezrKxTCZTErRzNXS8LhWu68DVPgsJZPJIJlMAoCm6JHnSKIxcW0W28vOTA+9tbU1RKPRrCos5Rgzk8kgk8nkNIe12uKzimoVWcD6YJNi3IBGezrKxTESiUj3YqU3PKbkhwqfBYiiiGQyidnZWZw4cSLvDUusrmJcm1rn1IPjOEQiEdhstpwqLIWgJ2CBQAALCwtoaGjI6n/n9XqRyWQQi8UgCEJeF1SxWJk6YSUHbY/PLOXe/1KKYyqVgsfjyeopaabhsVoah17bqmpdIFUKVPjKDM/zUtRmMBhEf39/3tfY7XZEIhEEAoG8rk2jnQzyCR8R2fr6evj9/pJMKmrCJ4oilpeXEYlEMDY2lvXZyCp7aWlJSvEIh8O6wQskcIH8W6uum/1GEISqEj4j45lpeKwmjul0WrOnIwkiK2dkZy3fK1T4ygQpKExcm0bbynAch3v37iGdTuPixYuGcuWMRJzpWV8bGxtYXV3FwMAANjY2DPf4yzf5KMfkOA6Tk5Ooq6vD8PBwzmcjq+yGhgbNnDBl8IKRhGkS8q5WasvKKjFWQi2+0lJocIt8oaaFWk/HaDSKdDqNjz/+uGwNj6nwUUqKKIpgGAbpdFp3L08Jsbo6OjqQSCQM3WgkTSEfRCDl8DyPubk58DyP++67Dw6Hw3T0p1HhSyQSmJycxLFjx9DZ2Wno3FrP5QteyFdqSxRFOJ1O+Hw+xONx2O12aa+xnB0MAOAX//dbWY9nvv25soxjNVT4Cketp2MikUAmk8Ho6CgAcw2PjfR0rMbFnhmo8JUYnueRSCSkhFcjk4EoitjY2MDa2hoGBgbgcDgwOztraDwz0ZpyMWEYBhMTEzh8+DB6enqk91mqLuzAp1bh9vY27t69m+O21YvqLGbfzUipLVKkeW5uDgBUOxhoTSI+n8/0CvvMC29p/v3D331AsxVRoVCL72CPpxRaMw2PjfZ0HBwchNvttuTzVBpU+EqEvFksANXVoZp7UB61SawulmUNiRkZx0yaAgDs7u5icXERZ8+ezZlwi21GK8dmsyEajYJlWYyNjRkWCyuiOkmR5rq6OnR0dOQ0h1WbREglEfkKW+l+krtVyTWgJXqEi//pOv7n88Cjjx5c66/WhMiK8cwursz0dDRScrCaocJXAohrM5PJZFVgUUJEhdywxLXp9/vR3d2dc5wRzAgVy7JYWFhAPB7XrMJSKosvk8lIVuv58+dNWQOVkMdndBIh7ieGYZBKpSSXKonq+8oPk4bexx/eBh59tKiPkgO1+EqH1cJXrgLV5Lp2Op1U+CiFI3dt6okesOeCI1FactemMmrTTMCF2t6dGqIoYm1tDYcPH9atwpIv6Vx+nJbwRaNRqS/gzs6O7neiNmEelMot+dxPDMMAP3zH8PnOvPAW3viNkzlh74VM8NWewE6FtnhocAvFNGS1PzMzg7a2NkN7NA6HA5lMBktLS7DZbJJrsxiMBLeEw2EsLy+jubkZfX19ec9HoiP10BInIuhDQ0NwOp3Y2trKey6j5y4H5bQuR/9X46JHeOutBTz88AmEQiEpQApQb+1D3KpqHdOrXYgAayfu/bD4aPf18kG/2QIQBEFybfI8b8hCIq+bmprCsWPHslybxaDn6hRFEaurq9je3sbJkycRjUbznq9QV6cgCJibmwPHcVKFGZZla7Y7w3/7b/r7elp8dxX4zRMncv4u717AMAzi8Th2d3c1O6YDexZnNBote6QqsD/CZyX7scdXy67IckOFzyQcxyGZTEpRm8R9mY/19XWEw2H09fWVTPQAbeHjOA7T09NwuVwYGxtDLBZDOBzOe75Ci1lPTEygo6MDfr/fVISo2mRZDRbfczdKez4j3QtIGgfDMAiFQmBZFnfv3pUiVUVRlCJV5UE4yqTpQqh24bP683EcV1S5QCNU8/9XPqjwGUQetSnvqJAvuITneczMzAAAjhw5YmpyMXKzkaAVOWqlzgpNe9AbVxAEhMNhzM7O4tSpU2htbc05Jp+oHOQ9vnJx5oW3Csrvk6dxuFwusCyLwcFB6Xm9cHeGYfImSvt8Pk0rpBZa6FSza7XWoMJnANIslmXZnNw8PYtPGbW5srJiOlozn5/f4XAglUpJjzc3N7GyspITNGPGkjNauWVnZwfJZBLDw8OqlSkKnSgOuqszX/qCVSi/f7ORqvJEaRK1qlV/kpTlcrvdVS+AVkCFr7xQ4csDKTsmiqJqQrrT6cwSHgIpA3bu3DkpkTpfZ3U5ZoRPEARpj00rZ66U+Xk8z2N3dxdOpxNjY2Mln+hq3eIDCrf6COWKVCXnltefJPvdCwsLUmUcec87pVu13N3SqwHafb28UOHTgNzcDMPAbrdrrr6UgiJ3bSqjNp1OpxSIkA8zrsl0Oo3x8XF0dnZm7bEVcr58okMqvtTV1aG1tbUsE5iV9TMrVfgqGbX6kzs7OxgeHpZc+cqaqmTfkTwG1CNV5eJYyxMztfjKCxU+FfRcm0qcTqckKKQepTIhnVCOxHQS3TcyMqK5QgeMW1F645JWQmfPnkU8Hi+JVaa2x2dlq6BSC1+luDn3O53BSE1VeaSqVsFxeaSqXBRFUbTsM+7HwqiQyi1GoQs9Knw5sCyLlZUVZDIZHDlyJO+NRdyXaq5NJXKRzEc+4RNFEUtLSwiFQmhpadEVPcC4mGi1Erp79y5CoRBGR0fhdrulVkv5yHeTVWtUZyko1t1pJYUEt5iJVJV3LWAYBgzD4L333suKVFUrG2e0K0qpP1uxlKtyC6GWu68DVPgk5K5N8ruRC8NmsyEcDsNms+WtR1nIHp8amUwGk5OTaGpqwsDAgFRoOd/7NPJ51FoJTU1NwefzYWRkRJoAzFivZrHa1Ukb0VbueGoFxzOZjNS2SxTFrDQOvcLMaikcRlr67Ifbkbo6ywsVPnzq2iSrLJfLZUigiGvTbrfj3LlzhqzDYi2+SCSCmZkZnDhxAu3t7aYKWhtBLjrk8x09ehRdXV05xxkVDLOT4kF1dVaKm3O/sEpo5RaYzWaTCo43NjaqHi8vzKzX0ketUzrZx9yPYJxatsjKTc0LXygUkiZ7UmvTSFI6cW2ePXsW09PThi5SM8KndIuSWpubm5s4f/48fD6f6XOaYWdnB8vLy5quWzPpEVrCx7Is1tfXVSP+rBS+SuYf3noLX/yceXdnNSeUm3U9ytM4jEaqErdqKpVCMpmUXKtKcSRW5EGLVK3Wa8MoNS18oijil3/5l/GTn/wkJ/pSS0x4nsfs7CwEQTDVaoec16irU55PR1oXORwOqRyY/LhSugXJBLCxsYHR0VHNhPtiG9YSa7K7uzsnqEEURSSTSdy6dUt136aYCiNqVOoeHwD87r8BXzwY23yWUQ5R1+uUHo1GcffuXQwODuaIYzgczrpu5edRXreVFKlaKe9jv6hp4SOo5eapCRSZrHt6etDd3W364jHjVnM6nUgmk3kjRc2it1pmWRYTExMA8rcSMip8am5L0g9wYGAALS0tqu/nnXfeQV9fX04oPMMw4DhOWsUr923IY6Or71ItHMbHx4s+RymhFl9px3M4HLDb7fD5fJK3RQ2e56VYAa1IVZLGobbnaEVj2Epe6FlFTQufXt88pfAZidosdDyt9xCNRrGzs1PUmMpzak0asVgMU1NTOHHihNQ9Qg+zrk5g74ZbXl5GJBLR7AdIsNlsaGxs1N23USZRE9cUqTAizxMjk4xaR4NSTASP/fdA0eegGGM/WgQZHc/hcORN41BGqsZiMamPIykAwDBMlsdDfu2WIlK1WhdFRqlp4QP2LCuWZbNWWvLJ2ohrU62zejEIgoD19XUkEglcvHixZPk8WmK1ubmJ1dVVDA0Noa6uDktLS3ktBrOuThId6vV6MTw8XPR3peeaIijzxJQdDeS5YIIgqE4wlUAhaQ3U4isdpY6wVItUlZNOp3Hz5k3J45FKpRCJRLC5uZkTqarVqqpSrt1Kpea/nebmZkQiEXR0dGT93WazGXZtKjurF0MqlcLk5CQaGxul7g9GMDLRKQNhBEHA/Pw8MpkMRkdHpbGIS7ZUwpdMJrG0tAS/3y8VzbaCfHlioihiY2MDwWAQzc3NSKVSCAaDWaHwWqW3fD7fgQtoqBb2w9VptdCSKNV8kapyj4f82lUWHJdfu263Gy0tLZZ9nkqECl9zM6LRaI7wkVw5I25GsidoJOhCzzoMBoOYn5/H6dOn4Xa7sbi4aOgzkHPmW5WqtRJqb2/HqVOnskTOiJAbFT6WZTE7O4vBwUFDzXqthITCu91udHZ2qh4jz+9UBjSkUimIooh/+BFj8Ts3BrX4SkMl9uIzWnBcGYwTi8WQyWRw3333leOtHxio8H1i8RGIa5PneYyMjBjabDYTrakmKvLKKCMjI/B4PEin06Zy/owIHxmb5AKqtRICjIlavmNIE9x4PI6zZ8/q3pz7Sb6AIyMu1Sd/VHn5e/v9vZaTahe+UhWo1rp2q7Fgg1lq3k8jF75EIoHx8XE0NzdruhjUMJufJxdJlmVx69YtcBwniZ7Zc5opQL25uYn5+XkMDw+rih45rhjh43keU1NTSCQS6OzsrOgKFJVcskyO2eR4usd3cMejVVvKD7X4PhE+0seOuDaDwSA4jiubxQfs5QdNT0+jr68vx9VaauHjeR6hUAgejwejo6O6N5aR82lNqul0Gnfu3MHhw4fR09OD+fn5il5hHhThA4D3338/Z59RnttYrUKnZD+ESM/iLzXlbklU63U6ASp8qKurw1//9V/Dbrfj8uXL0gVXqrqaasdyHIe1tTVsbGxIkZRKzKY+6I1PWgl5vV4cOXLE1F6gGYgL9fTp09LmeTF5clZYLQdJ+C5cuJBVk1Ie6ceyLABIuY0cx0EURQQCAUkcq8WKqHYLrNzjHZTrvZzUtPDNzMzgr/7qr/DAAw/gi1/8YtbNZMaKM3OszWbD0tISvF5vThWWQtHrmk43GRpiAAAgAElEQVQCZs6cOYNIJFKyZrRK1tfXsbGxkdONvVARJa8r94RTrPBZWZ8zXxi8vHv6+vo6ksmklB/GMExObmOpKotYPZFacV0ox7NaaGk6Qnmp2W83lUrha1/7Gp5++mns7OzkXNhmWwgZaTBLJqL29nacOXOmoPetNb5aK6F79+4hGAxKrYTi8bhh4TP62Unnd47jVF2ohQqfVZbYQbL4/ur/eQtf+4/a+Xzy7unJZBJerxd9fX1ZxygbxKqVi3M6nbouVSX70Qmimi0+juN0q8OUAurqrFG8Xi9+/OMf48c//jG+//3v5zxv1uJLJpO6x2xvb2N5eRmdnZ0lqcIiR+nqlCeLy1sJ2e12Q59JTUjVEAQBN2/eRFtbG3p7ewvqraclOlYVqj5IwvefbwFf+4/GjtX7XvNVFlHmh+UrF+d2u8HzvCS25RYl6uqkFEvNCh+wN+mRPD4lhQasKBEEAQsLC2AYBmNjY9jZ2TFsTRmtCCMfnyTd9/b25iSLOxwOqV6gHkZEJxaLgWEY9Pf3o62tLe9nMAsVvv3DSJsfeX5YIpEAy7KYmZkxXS6uEKo9qrOcwS30Wt+jooTvjTfewDPPPAOe5/HUU0/hueeey3r+7bffxrPPPotbt27h6tWrePTRRwEAN27cwNe//nVEo1E4HA784R/+IX7t137N0JgtLS1ZeXwEp9OJVCpl6BxaIkmSxNva2tDf3y+1PMpnHRKMVoQhrsmdnR0sLS1hYGCg6FZCeseRCNi6ujrNlAj5uUjghRmsakZbjPBVcv+9croflflhmUwGwWAQY2Nj0jFGysWRkltqblW9iV8QBEtdddVm8dGozgoSPp7n8fTTT+PNN9+E3+/HxYsXceXKFQwMDEjHHD16FC+//DK++93vZr22rq4Of/M3f4NTp05hfX0d999/Py5dumSoLE9LS4uqxWcmqlNN+EKhEObm5nKSxAtJU8hXEYbk59ntdt3iz0bH1jpOFMUs6/XWrVt5V8PF7PFRi+9goCayRsrFcRyXVRHHaLk4I/vppcTqYBrq6iw/FSN8169fR39/P06cOAEAeOyxx3Dt2rUs4Tt+/DiA3G7Ip0+fln4/cuQIOjs7sbOzY0j4mpqaEIvFcv5eqKuTBJUEAoGshHS1Y82cVwuWZbG6ugqHw4Hh4eGS1ddUfnaWZTE5OYmmpiYMDQ3BZrMZTnQvRFisdHUeJIwWrLYy4KSQYBN55/R8JbeU5eICgQBCoRBmZ2fhcDhUg3B8Pl/JchtLVYfXKOXO46NUkPCtra2ht7dXeuz3+/HBBx+YPs/169eRyWRw8uRJQ8c7nU7VCdZMVCcRSZZlMTU1BZ/Ph9HRUdWbpZT5gaSVUEdHB3ieN12kWgul6JB9w+PHj2cl2peitJkWVlpi1OIrjnK5HrVKbt25cwe9vb1obm6WXKpEHGOxGLa3tyWXKvBpoWY1cTRiWVVbMM1BW+yVg4oRPrXJx+x/0MbGBr785S/jlVdeMX2hKlfIZrulcxyH8fFxHD9+XLPgMTlvKSw+ss82ODgIjuOwublp6H0aESH5uHr7hqWq6VnMey0Wq/YSrcZqi8/KyVQuREZcqiRKlYhjIBCQficLRr3cRsBaseA4jgpfmakY4fP7/VhZWZEer66u4siRI4ZfH41G8YUvfAHf+c538OCDDxp+nbwZaaHCt76+jlQqhQceeEA3TBwo3tVJWgml02mpP2AsFitq704JCZZZWlrSbRpbrPAxDIPJyUnVKECGYRCLxeDz+YqOAsxHIcJXyYEtVrOfwpcPI10MBEGQolQZhpFSOBiGQTqdRjwet7xcXLkszGqu4WqGihG+ixcvYm5uDktLS+jp6cHVq1fx6quvGnptJpPBI488gieeeAK/+qu/anpsMtHKV41GJnXSyUEUxby5UYRihC+TyeDOnTtoa2vLaiVkNO/OaGK6KIrY3d1FZ2cnRkZGNG8UI9+RVpAKKW929uxZtLe3S4nVZCXOcRy2trawubkJlmUhimJObzGfzyf9XuhEcRCDWxYXF6W9cC1qxeIrBXa7XbqW1KKU3333XVPl4tTEkQarVBYVI3xOpxMvvfQSLl26BJ7n8eSTT2JwcBDPP/88Lly4gCtXruDDDz/EI488glAohB/84Ad44YUXMDExgb//+7/H22+/jUAggJdffhkA8PLLL2N0dNTQ2KRQtZa7RA1S/7K7uxtHjhzBRx99ZGgCMOPCkwsfEYr+/n4cOnRI87h858s3NsMwmJubg8fjQX9/v+6xhVp8GxsbWFtbw/DwsFShQplYnUgk0NraKu0pEpcVmXgYhsHOzo70mEyGypwx8q/WqvwgCt9/eGURM9/WFz4rsbqSitV7boC5cnHkmrSiXBylMCpG+ADg8uXLuHz5ctbfXnzxRen3ixcvYnV1Ned1jz/+OB5//PGCxyVJ7EZdq2TfS95nzmhtSbPFp1mWlQpanz9/XrWUkVHhyzfJk7qefX192N7ezns+s1GdJB0ilUrlrVOq3HuTu6yam5tVXyPPHWMYBtFoNCfQQRke73K5wHHcvkym5cRKMbI6r87Kz2Z0USQvF6d1fRotF8cwDGZnZw2ViysEKrAVJnz7RVNTk2oSu7JyiiAIWFpaQjwez9n3InuCpXRp2Gw2bG1tob6+XlcojLowtS54URSxsrKC3d1djI6OQhCEkgXLkGM4jsPExASampowODiY9+YrJI/PSKBDJpPJshrD4TCSySQ++OADaeGidKPKRZKSy0F3dVo1lpFycel0Gh999BFaW1sNlYtTWo9G3isVPip8APSrt5Acnkwmg4mJCbS0tKjmyxHhU+btFUoqlcLy8jK8Xi/OnTuXNz+vUHiex8zMDBwOh5SCkclkSlLhhRzDsizGx8dx7NixnIhXK6M65atyAs/ziMViUkAUx3GqeznyiecrPzRWeWc/oXt8pRvLyv05ImzK/pwEZbk44tkgj/OVi6P5gXvQbwHIW68zkUhgdnZWdX+NYCZoJV8NTlL1xe/3I5lMlm1SSaVSmJiYkJrGEgrN91MjGo0iFovh/vvvN9XVfr+6MxjZy8EP/7+yv698XH7hLfyTgUR2K9iP7gxWjWd18nq+lkRauY3Kc2iVi3M4HLh48WI53vqBggof9oQvHA7n/N1ut2N1dRXRaDSnz5ySQiq9KG8ouctxZGQEmUxGtapMKQiHw5idnc1qGksoVU1PsjdZV1dnSvTMvIdiMSuwleImWsjzvJXiYPUen5XsR0uiYsfTc/nTdIY9qmc3vwhaW1tzBIbjOMmFMDY2pit6QPH5eWQPLJVKYXR0FB6Px9Q5zbC2toaFhQUMDw+rlnUzemNoiZMgCJidnUU4HMbo6Gje86kJDy1ZdnCwOqrTSqqtastBi2AuF9Tiw6fpDIR4PI6pqSk0NDSgu7vb0IVv1uKTH5tMJjExMQG/34/u7u6s48wIX77VHEnUDYVCqk1jzaJV03NiYgKtra1SrmGhtTqNfp+1yvj4uGbbH8A6Ua9mK2I/OjPQfbjyQ79hZHdo2NjYwOrqKgYGBhAIBEpWV1OOvGzZ7u4uFhcXce7cuRx3oJlz5kunIMnvTqcTp06dKsnNrHx/WjU9C8Gq7gxmqaSKLWfOnJGCcBKJBHZ3d8EwDFiWRTqdhsvlQiAQyAp0KDbhXw0qfKWj3OXKAOrlAKjwAfjU4nv77bfR1tYmlQKLRCIl7cJOIBbf4uIiotGoZkmwQtynajcNKWbd39+PtbW1kgmK3B0ZCASwsLCg2QtQj0I6t1OAX/jTdzGuEeAyPz8Pr9eLxsZG1YR/IlZmE/7VqGbhqzZXJ2UPKnzYE4br16/j537u5/Cbv/mb0k1cSNNYoywuLuLQoUO6JcHMJrurjU+KWQ8NDaGurg4bGxslFT6e57GysoKdnR2Mjo5KbjYzqE2cB7GiitXoXZmk6XFzc3PBCf8ktF5pLZIfIgjVLHz7YfHR7uvlp+aF75/+6Z/w3HPPobm5Gb/1W7+V9VyxzWjViMfj2NraQnt7u+HWSUZQCp+yaSy5mcwItJEJLRgMorm5WbMNU6FQi6/8FJLwH4lEpOLNxLUuiiIcDkdW8n+5Ev6tnrz3Y4+vVLnAatDu63vUvPBFIhG8+eabuHTpUs5zpWohRNja2sK9e/dw5MiRkk8KcqFQaxqrdpyR8+ntGc7Pz8PpdOLs2bO6N1MhNxoVvuIohRWmlvCvhOM4LC8vSxaiPOGf3A8ejyeneLPP5yuoPmW11wXlOM5UzWBKYdS88H3pS1/KqiVZaGsivWMFQciqUUnq8xnFyCRGhJdEpGoFmJhNTlcTvng8jsnJSfT09CAcDhuavMxOxJXo6qykwBYjWLGydzqdcLvd8Hg8WY2kCcou6qRMHLEaSX1Krb1GpdtvP/bcrCxVR/f4rKHihO+NN97AM888A57n8dRTT+G5557Lev7tt9/Gs88+i1u3buHq1at49NFHpedeeeUVfOc73wEA/NEf/RF+/dd/3dCYWhNEKYSPlDprbW1Ff3+/tPdipsmtkbJJDocDwWAQwWAQAwMDmqtGM81o1Y4jBboHBwdht9sRDAbznouImF6tUOVz1OIzxpkX3sKMSoCLlYsGPSvMSKURjuOyhDEYDEp7j6RMHDmH2+0Gy7IIh8OWdDXgeT5vDm+px6NNaMtPRQkfz/N4+umn8eabb8Lv9+PixYu4cuUKBgYGpGOOHj2Kl19+Gd/97nezXhsMBvHtb38bH330EWw2G+6//35cuXJFtb+WGmoCU0g1FjnRaBTT09M4efIk2tradI81876UiKKISCQClmU1I0TNjq0sfC2KIu7du4dQKCSNQfZ5jH4GvclR6zWUyqdYt6rT6URjY6NmdR+Sf0qaEwuCgLW1NTAMg0wmA1EU4XK5NK3GYoRkP1yd5czjo8K3R0UJ3/Xr19Hf3y812Xzsscdw7dq1LOE7fvw4gNzCzD/84Q/x+c9/Xqql+fnPfx5vvPEGvvSlLxkau6GhAbFYLKuSidGuB0DuBbW+vo719XXVVkJmc/44jtMUM47jMDU1BQDo6enJ65YppA4nz/OYnp6Gy+XC8PCw9N0bFadCcvIq0dV5kLC6SHU5xUHeKNbr9SIYDGJwcDBrfKXVGAgEpMekPKBa4Wafz6eburEfwS3lGq+ao2/NUlHCt7a2lrVP4Pf78cEHHxT82rW1NcNjk1w+ufAVcpGQcl2CIGi2Eiq2vBmBVHzp7e2FIAiGrS8jVixxdabTady5cyenkDU5V7E1PYPBIO7evSsFQJBJyWaz0cotBwQra3WqWWA2mw0ulwsul0vqj6n2OnnXDWXCPwBVqzGdTlsqFrRyizVU1DestsI3k0Rb6GuBvZ58ah0azCAIAsbHx9HV1YWenp6S7B1qCR9pGksqvmxtbUk3cL7zpdPpvMfZ7XbEYjFMT0/j1KlTqi7jYoVvbW0NW1tbOH/+PERRlCalcDiMeDyOcDiMd999FwBUV+r5Vuul5KAFtlRrW6JCXY/5euGJogiWZbOsxlgshlAohFgsBrvdnpPwr2wSW4rvwIrKLZQKEz6/34+VlRXp8erqquGu6H6/H//6r/+a9drPfOYzhsfW6tCQLzCDEAqFwDAMhoeHNVsXEYqx+JRNY0nCeClbCQEAwzBSlwitycJMMWv5wkQURczPzyOTyWBkZAQejwd2uz2r4ks6ncbt27dx4cKFrD0ehmFUV+vKRGv5pFTt7h2tABerOAjClw+SrO92u7MS/lOpFE6dOoX6+vqchP9IJIKtrS2kUinN65D8Tq7xfJTbbVzt94JRKkr4Ll68iLm5OSwtLaGnpwdXr17Fq6++aui1ly5dwh/8wR8gFAoBAH70ox/hT/7kTwyPLa/XKYcIipb7QRRFrK6uYmdnR3eDXo6ZwA25oJG9NqfTmZMwblT48h0niiKWl5eRSCTQ19en2y3aKPLPy3EcJicn0djYKEW55nuNfI9H6z2rdVYnIfPA3udWiiL5txqhFl9pkO+5FZvwT0RNLQjH5/NZ4uKkwrdHRQmf0+nESy+9hEuXLoHneTz55JMYHBzE888/jwsXLuDKlSv48MMP8cgjjyAUCuEHP/gBXnjhBUxMTODQoUP44z/+Y6nJ4vPPP5/X8pKj1YyWVG9RuyiVQnTnzh1DeT+FlCLTahpLKMV+G8/zmJqagsfjMWxpG4GMmUqlcOfOHfT29qKrq0v3NWaCW4wkWrMsK01GxJolj+PxON57770cUbRyQjrIWN0RvVIT2I0m/Mv3GpUJ//F4HB9//HHJEv4p6lTcHX358mVcvnw5628vvvii9PvFixexurqq+tonn3wSTz75ZEHjKlsTEbT245LJpJTETVoJmdm7M4rD4UAsFsO9e/dUm8bKjyvG4iOidOTIERw5cgSrq6slrekZi8UwMzODM2fOaNaOVL6mlOkMJPhBzSJ/7733cOHCBWm1TlbqygokBwlq8ZWGUkdZOp1ONDQ0aBZyf+eddzA4OGjIe2Ek4Z+iDv2WPqG1tVU1ClStbBnpRHD27NmsKDIzJc6MEo1GEQqFcN999+m65YrZ44tEIpiZmckS1lIKD7GwRkdHNd2VSqxMZ7DZbLDb7ZpurIMW2GI11Sx8+2HNFpPwT+YAtVZUJOGfQoVPwojFR/a/IpGIaicCM0WtAf2biqRFMAyDrq6uknWAVx63ubmJ1dXVnHxDu91uKEpUD1EUcffuXWm/0KjoAdbuRVRDzqAywIVafKXDqs9m1Lo0mvBPhDGRSCAQCMDtdmflRNcyVPg+QWuPjwgfCcqoq6vLSuKWU0gzWrXzkKaxHR0d6OjokAJ29DCabE8sOVEUsbi4iGQyidHR0RwXiVbJMiUkOV35OQRBkPY/Dx8+rDtZ7bfoVIPw7SfVLnxWUaqqLfJgMHkakpX5lpVOdV5BBaAX1ZlIJDA+Po7Dhw+jv79f88YrRX5eNBrFjRs3cPz4cfT29hp2nxoVKmKV3r59GzabDUNDQ6o3WzHBMplMBjdu3EBTUxNOnz5t+L3tF9UqfNUoRtUsfLRAtXVQi+8TtIQvmUxie3sbY2NjeduFOJ1OZDIZQ+OpCZ+a29GoFWl0kiNVK3p7e3H48GHN4woVvkQigYmJiaz6pJVed7NahO+dd96B2+2Gz+dDNBqFy+WCIAg5jWNLjdUWn5XdEqyECp91UOH7BKXwCYKAxcVFRKNRHD582FCPrEJcnWQs0rZI6XY029ldj3A4jNnZWXi9Xl3RMzOuXNRINZmBgYGsqLWDLHwHKbCFtKJiGAYzMzPgOA5bW1tSIASArGaxytSNQifdat1PtHoxRLuvWwcVvk/wer2StUZaCbW0tODkyZPY2dkxdA6zrk6O48CyLCYmJtDc3JzTNJYcVwrhW19fx8bGBkZGRnDr1q28x5u1+NbW1rC5uaka9FNIkWorqZZ9j8/+lynMfLsHbrcbdXV1OHLkSE7qiFpEIHlM9oC0ysM5nU7V76pa9/isLupsRUuiarnWi4UKnwKSqkBcdfF4vOh9O61jE4kE5ufn0dfXp9o01uw51SDWZDqdxujoqOEby0znheXlZQiCoHn+g2DxVfL7KyVGIgLlif7RaBTb29tgGEa6D+Qd1X0+HzKZDDKZDHw+X9knViuFz2rXY7lbElE+hX7Ln0Bq9f3Gb/wGXnvtNalUl5kUBTMWH0lMHRkZ0XWjmhUN+SqV4zhMTEygqalJtzyYGkYEl+M4RCIRtLa2YnBwULf8mJ6rZb/dMNW4Ci7UWjFSzFleN5X8zM7OguM4iKKoWa+yFPuM1Sx8dI/POqjwYc+1+eyzz+LevXt4//33s256s13Y84mFKIpYWlpCLBZDT09P3r3DQsTK6XRKLYuOHTuGzs5O1fehd+58gksqvdTV1aG7u7uoc2lhtEB4KVATX60KQQeFcnxv8m7oJFQ+EAhgZGQEbrdbtctBNBqVEqzl9SrVOh3ks3iqOYKU5/mcbYJSUo0LvEKhwgfgH//xH3HmzBmcO3cuZ4I2U40ln3UozwU8duyYofZAZiDCF41Gs1oWaR2nN8noiRXpLH/mzBns7OzkFbVihC9f9/lSoGWRPvx/TJd13HJw4YW38NG3P2epFS1fnGh1OZAjr1cprzzCMEzWPqNatw2t3NdysB+uTjNFHsxChe9TqPABeOSRRwAAP/nJTxCNRrP228xcLHoTPLHAjh49iq6uLuzu7iKZTBb3xlXGX1tbQzgcVg0yMfI+8x2zvb2Nu3fvSikXgUCgbMKXz0VaSvbb3VoqYvswplmrPF+9SuU+YywWk/YZSY+8uro6zfY/pZrgq8nVSbuvZ0OFT4ZW2TIzXQLUIAEzcgvMbHmzfBeuIAiIx+MQRTGnZZESI/t3yrFI+bFwOIyxsTHJWixGRPNhVdBJteTxyTnIKQZ6+4w3b97EiRMn4HK5VOtVkvY/pJu6mjvVqMVopXVJxqPBLdZQsd/yG2+8gWeeeQY8z+Opp57Cc889l/V8Op3GE088gY8//hhtbW147bXXcPz4cbAsi6eeego/+9nPwHEcnnjiCfz+7/++oTGbmppUhQ8o7OYmTWMDgUCOBWbGhUqEQ2s1yLIs7ty5A5fLhaNHj+a9Wc0KkSAImJmZgd1uzynXZuRcegLGsiymp6elqED5jyAIYFm2rPse5P1Vm/BZidXpDCQXUat+LdlnlFuN29vb0j6j/BxqPRqJ+FjhZpdDu69bR0UKH8/zePrpp/Hmm2/C7/fj4sWLuHLlSlaB1e9973tobW3F/Pw8rl69im9+85t47bXX8Prrr0vdu5PJJAYGBvClL30Jx48fzzuuVvWWfMKjRBTFrFqVIyMjOWJUSBd2tfFJpZQTJ04gHA4bLltmdGwiqu3t7fD7/TkTnFGLT01YGIbBnTt30NfXh66urqyJinS0v3nzJoBPgyqU4liKfnlqwneQEteVkCCSg2rx6WEk4ES+zyjvniJHuc9IrjfShorkvJF0Hfm153K5yvJ5rcjjo+xRkcJ3/fp19Pf348SJEwCAxx57DNeuXcsSvmvXruFb3/oWAODRRx/FN77xDekGTCQSUqKu3sWvJF+haiMXpd1uRzKZxNTUlNTbTo1ChE8JcaGSSimxWKzg1kRqCIKAGzduoK+vD+3t7QWfS68V0tmzZ9Hc3Kza3TqZTOLUqVOor6+HIAg5na1JvzyO42Cz2VSTrslEpUe1WXwjf/JTvPaI+v9Xuagk4TOCkX3Ge/fuIZFIwOVyIR6PY3d3FwzDSF1LSHk45TVX6D5jufP4qPB9SkUK39raGnp7e6XHfr8fH3zwgeYxTqcTzc3NCAQCePTRR3Ht2jV0d3cjmUziL/7iLwx3Ym9pacHm5mbO3824JQVBwJ07d6QJXYtiClpruVCNdmgwIrpkBXzhwgXNyYGMme9zKIVPGSBj5HV6/fIA/cRrlmV1K5IA1RPcQqjWYAarUgzsdjscDgcaGxvR09OT87woilmNi1OpVFbDWLLPqOVOVfsMNI/POipS+NQmIbVgC7Vjrl+/DofDgfX1dYRCIfziL/4iPve5z0nWox7Nzc2Ym5vL+btRkVpbW0MqlcL58+fzdhk3s88mFyqy32az2XJcqEa7IOQbm5Q38/l8efMMjYiovBWSXLDzWWJmojrzJV6THmVkoorH49jZ2UEqlUIsFsPu7i7q6+uzxJBSeVidwK51jdpsNng8Hng8Hql5sxJlPuPOzo70mOQzyoUxnU4jkUigrq6OBrmUmYr8dv1+P1ZWVqTHq6urOS5Dcozf75cqiBw6dAivvvoqfuVXfgUulwudnZ34hV/4BXz00UeGhE9rjy+f8JGmsYIgoK2tzdCqrZDEdNKnr7OzEz09Pap1PY0ItJZYiaKYVSz75s2bReX7EWw2G3iel74jtT1PrdeVKqpTnjStZGZmBocOHUJjYyMYhkEgECjJmPtJtVp8lbafqIfL5YLL5dLcauF5PksYeZ7H8vKy9DsATS9FIfuM1Xg9FEpFCt/FixcxNzeHpaUl9PT04OrVq3j11Vezjrly5QpeeeUVPPTQQ/j+97+Pz372s7DZbDh69Cj+5V/+BY8//jiSySTef/99PPvss4bG1Upn0BOUdDqNiYkJdHR0wO/3Y35+3lSaghEcDgfi8Tjm5+fR39+v6bp1OByGkuLVhI/neUxOTqK+vl4qP1aqVAVBEJBMJtHe3o7jx48bvgGtqvFJ3g+JFHzwf/u47GOWmy9dC+LGuf1+F+Wh0jqiF4rD4cjaZ1xZWcHw8LD0vNJLkUgkVPcZ1dypXq+XCp0OFSl8TqcTL730Ei5dugSe5/Hkk09icHAQzz//PC5cuIArV67gq1/9Kr785S9LQnD16lUAwNNPP42vfOUrGBoagiiK+MpXvpJ1MenR2tpqyuIjFUzkYlTKNkIEUtdzdHRU05UHmOvCLv88JAq2p6cH3d3dWccVK3zk3E6nE319farHaLkzrUpgr7bgFgKd+Ipjv/fc9LwUQPY+I7Ea5fuMwN585PP5pEUnZY+KFD4AuHz5Mi5fvpz1txdffFH63ev14vXXX895XUNDg+rfjaBl8akJ38bGBtbW1nICNMwErQD6rhuSNB6LxeD3+3VFDzDXhZ3cGLFYDFNTUzh9+nTOXoWZ/Ts14vE4JicncerUKdW9U2BvVUvCx+12e1brFJrATtlPrG6BZBb5PqMWJJ+RXt/ZVKzw7QfNzc2Ix+M5f3c6nVIjT72mseRYs9GaantoPM9jenpaSkonrg0j58sHEaudnR0sLy9jaGhIVVSLsfhIqsXg4CDq6+tVxV0QBCkVgeQ+EnieRzgcRmtrq3QMkNtTrBRWTbUKH7X4isNKi69cVVvIPmOttN0yinX1eA4AWmkLRMxYlsWtW7fgcrkwNDSkeqGa7cKuJpLpdBo3btxAS0sLTp8+DZfLVVTQihK73Y5QKIS1tTVd92mhwre+vqTnY1UAACAASURBVI7l5WWMjo7qph9wHAe73Q632w2HwwG73S65YW/cuIGuri6pA4AgCJJ1yHGc9MOyLDiOA8/z4Hleih6V/+RDLhAHOXFdyfD/8m+Ix+Mld73XClaWLCtn1RZyD9CF0KdQi08FpfvR6XSCYRiMj4/jxIkTmsnc5NhiGteSfcNTp05Jk34p8/MEQcDKygpYlsXY2JjujW1kXLnwiaKIxcVFJJPJnKa0cgEioudwOOB0OrPeQzKZxK1bt9Df35/1PWt1/pafUwul+1R5rmpuRDs/Pw+GYaRrWq3yjc/no/ljKlhZsoxWbbEWKnwyyMWhFL5IJIJwOIz777+/JHltWsdubW3h3r17qvuGRvPz9MYm5ccaGxtzBEfr/Rm1+Hiex9TUFLxeL4aGhjRvNGKZOZ1OycojRCIRTE5OYnBwMCsEnBdEJDIcXA47fK5PJ4d8Lk9RFMHyAhIZDh6HDV6X+sRCKvxUo2U0Ojoq/a4Mnw+FQlhfX89pB1SOsnAHEStdnbT7urXQb1qBz+dDKpVCXV2d1DQ2EonoVg2RY7ZxLelavbS0hHg8ntX5gGAmPy9fW6S+vj7U1dVhcXEx7/mMujo5jsPNmzfR2dkJv9+veSz5rC6XK+czbm9vY2lpCaOjo1min0hzeHN6B6FEBoAN9x9rwXCPsRJ0wSSLN6d2kMzwsNts+MX+QzjZ8en/IfneGYbBiRMnqtbqIyjD55XkKwtHEq6V3Q6s2h8l4mwVVu/xUYvPOqjwKWhqakI4HIbb7Zby2oaHhzE+Pm7o9WYtPmKF+Xw+nD9/XvUCNRu0oiQUCmFubk5qi5ROp0tW0zOVSiGRSGBoaEjXBWyz2RCPx9HU1JQjeisrK9je3sZ9992XUynjh5PbWA4kcbjZiyavE9eXQ+hq9KCrSTuSDdgTtX+Z3oXdBnQ3e5DhBPzbfAAdjW40eV0QRREzMzMQBEHqODHw4o/zfCMHjzMvvIWZb3/O0LGFlIXb2tpCMpnEO++8k2Mxlrqws9Ud0a2uEkMtPuug37SClpYW3Lt3D2tra1LTWEB/D0mOGYtPEAQsLi7i+PHjWflzSowKn9rEQsqPjYyMSGHPRhPD841LCk17vV5N0SNu0MOHD2N1dTWnXBN5fOLECUQSKdicApp8LrgcdkxtxvDW9C58LjvWwikcb6tDvceJWJpDF/Y+SyLNIZTMQBCBjkaP5AqNpFjcXIvAYbehweOA027DajiFZq8Dp7sakN5aQlNjA06cOAGGFWCzVV9UZ6lRKwvH8zySySQefPBB3bJwmUxGCr9Xc6UaEUarhQ+wzlIqd0siavFlQ4VPQSaTwde//nW88cYbkuiZwahIETdSV1eXruiZOaccUn6MYZicQBMzFqSWiG9tbUmVJm7fvq16jDxyU9knkGVZ3L59Gx6PB62trfhgfgvv3o1hJ5aBHcAv+F1YYRyogw1sxo46jxPj94LoavJh1L/n6tyIpPA/bm1gdiuJMJNBd5MXXzjfiSavE6+8v4o7a1E0+xwIJHnE0xxYXsC7C0G0ugQMdDfgd/9DD96c3sXSbtLEN0uRI98PN5JwbaQSiZowut1uacFUjdDu69ZChe8TRFHEn//5n+NnP/sZXnjhhYJEDzBWbYQkvx87dsywJWdmH4WUH6urq1MNNDF6PjXLUBRF3Lt3T6oko+WeUUZuwmZDIsPB7bADAo9bt26hs7MTvb292ImnsbjAIuMCGloAVhDw04CA7iYnBtpcuL0ew9RGHGGGxfpuDNcXtnHIC3Q2uOFxORHMABkeeH8pgZtrEXiddqxF0hAFEevRFHgBcNpt8DntSLM8goIN0wEOf/rDefS2enH0kA//9R/VE+yrgT/6y7fwnf/JmLvTLGasMOIK9Xq9UsSyHGXHA9JZnWEYZDIZCIIg1atVCmOhrYAqBY7j8hZtp5QOKnyfsLm5iUAggCeffFL1BiIh78WsOOVW2NjYGEKhkGqJNLWxjSIIAsbHx3PKjxVyPqVlSIpxA8D58+c1vwul6K1H0viv795DMMGixWvHkC+EfzdyCp2dnQCAnVgGs1sJbMXSaKtz4VCdC2FewPxuCnM7KTjsIlxuN1yCDZEUDybDY4sBpoMZNHpY2G028IIAEUA4ziItAESuibzzvIg0z8NhAzhBxHo0jZUQg3o3wAs2xAx9IweT13eB75Tp3KW0JvJ1PIjH45idnUVvb68kjMpWQFo98iq9diXP85od5SmlhwrfJ3R3d+NP//RP8Zd/+ZeaZcuKSWjlOA4TExNobGyUrDAzLkwjFlosFgPDMBgeHjbcg1APucVH3n9LSwuOHj1qKF3B5XIhkeHwn95aQJhhUecEJlbCYLra8HDD3sSW4QT8dCGA1TCDnVgGm9EUvE47Gt1O7CYzEEUbMpwAThTgdTmQ4XjYbIDLtmflRdMiPE4RnAA47DZkBBF6u5e8CPAAWHbvqHAK+FQeq5drN1YRTfE40V6PC8da4Haqp4WYxepO76T3plrbL1EUpVZA8gAceSsgl8ulKYzye9vqSj50j89aqPApaGlpwdLSUs7fSdCKUXeEfEJgGAZ37tzJCpaRn7MUkPJjDQ0NaGxsLMk5ifClUincuXMHvb29ui5gUjlFnq5wN8BgYScBUeCRzqTR2dKIXYbH7FYCa5EAtqIprIdTaPI6EU5y4EUBKVZAOJlCk8+FzkYPIIqY20kgkuQgALBhr+SQ+MkP88lXyArVL2CFcmtxA0c62/HeQgDJNIvPnM4tDqBM9Fc+r0YltQmy2Wxwu91wu92a/TCVwri9vY1UKiUJIynq7PV6wbIsAoGAbvPYUlHuqE4qfNlUtPC98cYbeOaZZ8DzPJ566ik899xzWc+n02k88cQT+Pjjj9HW1obXXntNqkB+69YtfO1rX0M0GoXdbseHH35oyJVQaE8+OfIanCSV4OzZszl9ucxYfGRfTq0hr7y56+TkZMny0RwOB1KpFG7duqVaxFr+3kjUnsvlylq5LuzGEUul4QIPn8+H9UgKbqcD/3RnCz2tPjAZAbPbCXicdtx3tBm7iQxSGR5bsRQYlsO9IA+H3QaWF7Lcl9WXal5e/q+bEbz4xcOoczswsZnAZ850wK7Y69W6bvSE0cqyXqWI6szXI4/jOCkidXNzU7V5rFrwjdfrLcpio+kM1lKx3zTP83j66afx5ptvwu/34+LFi7hy5QoGBgakY773ve+htbUV8/PzuHr1Kr75zW/itddeA8dxePzxx/G3f/u3GBkZQSAQMGypNTc3Fy185NitrS1sbm5mpRLIKaTKi/zmkO+5keauZt2neivBaDSKcDiMCxcu6HY2B4Dd3d2cHD1RFLGxsY06h4jdlAPbyRQEiGirE7EUSOJ4ex0ON3nQ0eDGUjCJOrcDHocdDY1O3A0kEcsc1IRyEXt2aeXB8SJcdhvsBi27fMIYjUZhs9mkQuJmLUYzWJHO4HQ60djYCJfLhfr6epw9ezbreY7jsnIZA4GA9LtaV3WjZeGoq9NaKlb4rl+/jv7+fqlz+mOPPYZr165lCd+1a9fwrW99CwDw6KOP4hvf+AZEUcSPfvQjDA8PY2RkBADQ1tZmeNxCmtEqsdvtWFhYAICcVIJCz6kUPpL43tbWht7eXunCNip8ep0hgL2u95ubm2hubtYVPY7jcPz4ccxvhPDuxxuIpVj01osYaHOA41gEozbA5kCjRwDHAw4RiKcFrIYYLOzEMdzTgnOHG9HZ4EY0w+NeIIl4mj3AogdUqujNrG6jrr4R/36gw/Br9AQsGAxK1XaAwixGtfNqYWUen9ZYTqdTt/qNXlk4Yh2rCSPLsrReqoVUrPCtra2ht7dXeuz3+/HBBx9oHkM2vQOBAGZnZ2Gz2XDp0iXs7Ozgsccew+/93u8ZGtdsM1olLMsiEomgvb0dp0+f1r2pjfbPI8cSQSPlx44fP46Ojg7N4/TQSmInkaepVAqDg4OSRalEHrnpqGvGB7tx1Ld24JDTjq1oClPhEHhHE7ZYAUkugWiKhwiA4wGOF5BgM/gfN7awvBlEo8+Nf3eiFS0NDfjbEIPdeP4u8hTz/N2NEH7w9VM41qaeZ2eGYDCIubk5jI2NSd6MQguJGxVGqyupFCJERsvCEXEk+bzRaBQfffSRalk4eZI/pTRUrPCpRVWp7W+pHcNxHH7605/iww8/RF1dHR5++GHcf//9ePjhh/OO29LSohvVqQcRpMbGRnR0dORdyZpxPxBBU5Yf0zrO6PnkyPP/BgcHpQhNJcp0hc3dOHhRRLPPBZZjwUZ3cJv3IcqKEAG4nQ6IIo8U/+n/lw17+3Q3t1icbrfh//7ZJtZiHFgOKE24D0WNv/94DQ67DQ+f7cD9R7P3bOe347izHsVGNA2fywF/qw8P9bWi3pM9TQQCAczPz2N0dFS3CSpQvCtV/jqWZaW0onK4UuWUS2S1ysK9++67eOihhyCKompZuFQqJX1+tYjUfNVvqKszm4oVPr/fj5WVFenx6uoqjhw5onqM3+8Hx3GIRCI4dOgQ/H4/fumXfkkqo3X58mX87Gc/MyR89fX1SCZzK3k4nU7VvxOCwSDm5+cxMDCA3d3dkkVrEhwOB7a3txEKhTT3DAHj5ciUx2UyGdy+fRvd3d3S96x2LhK5SdIVAMBpt0P4pCrH0so6BF8zFjeSqPc44XPaEUhkkOazFynkkSgCmwkBaVYAy1PRKzd/e30NDx1vwfxOAt/89/042bFnmUxtxvD/jm9gJ5bGSpiB1+nAma4GLO4m8RsP9sLt3BMBInpjY2Nwu91Fvx+jwphMJqWmycrFmNJKLEWzYisLVBPI+1aWhZOjLAsXi8WkyFRS/UZZFq6hoYHmCCqoWOG7ePEi5ubmsLS0hJ6eHly9ehWvvvpq1jFXrlzBK6+8goceegjf//738dnPflZycf7Zn/0Zkskk3G43fvKTn+C3f/u3DY1LVnlqPfnUxEwURayurmJnZwejo6Nwu90IhUIlbXEjiqIURDA2NqZ7QxZi8SUSCUxMTODkyZNZ+6HKCi9y0ZPvDZ5or0OLS8A7U2sICT5wEQbJDI9wkkMkxYLT0WEBQIThQDMRrGNxOwxOAP7LW3fw9YcOQ3S48fp4EBxsCCQyONzkRSzFwWYDAokM1iN7dVJLLXpGsNlsYBgGt2/fxuDgYI6Xw6wrlfybTxj3Q/iMUEhZOAAlyeutJipW+JxOJ1566SVcunQJPM/jySefxODgIJ5//nlcuHABV65cwVe/+lV8+ctfRn9/Pw4dOoSrV68C2Nun+53f+R1cvHgRNpsNly9fxhe+8AVD42qtDtWEj0RViqKI0dFRSTTN5ufpRVcS96PNZoPf7897M5rd4yOu04GBgZx9Cfl7IqtJl8uFpWAKP51fhyACP3e8BYdsCQx4Q1hsb8cR115k5j9PbCLEGNu/pKJnLVufOC7+eZ7BcnQNS8GMlAvpsgF9zYAAO0Iudq8kXMALe2IX6+vrkuilWB5up12KDi0XqVQKN2/exLlz51RTEApxpSoXdGrCSAJRDlqdS2VZOFEUD9xnsIKKFT5gz0V5+fLlrL+9+OKL0u9erxevv/666msff/xxPP744wWNS/YRlIWd5WKWyWQwMTGB9vZ2+P3+HOswlUoZGksvujKdTuPOnTvo7u6WrC0j58tkMoaO293dRTgc1nWdkmoYwF4B4ZVwCn93fRWNHgdsNhv+z3+dQYuDha+1AyvhKH6urxWrIQZp6rOseDgRuL21d62QqzclAtMhoN5tg9Mt4kKPC1w8iP9+K4io4MJP165jhxGR5B2o97pw+UwLzh5pkfaaShmAkU6ncePGDZw9e1YzIT0fWsJIfteyGDc2NtDR0SGlachfUwpXqpxyV4mhopdLRQvfftHY2Ih4PJ51s8mDW+LxOCYnJ3Ncg4Ri8/OAvfJjU1NTOHXqFFpbW7G+vm7Ykst3nCiKiEQi4HletfEtQRAECIKAQCCAhoYGOBwOTKzH4HHa0eR1YmdnBzuxDLZdHjzU5YTLbsPbswHUexyw6RYOO8hUbo5eMYiK3zmBR5IVEWR4XJuKg3G3o6vJgw9Wo/A4bHjoWDNiqRT+eT6BZrcNDmFvn4ksDrX68hlN0iaid+bMGc3CCaVAKWCiKGJ6ehp1dXXo7u4umytVTqW6VasZKnwqkGa0SuHjOA67u7tYWlrC4OCgZsNOIxGgyvPKLS4yxtDQkLTJTZrW5iNfioQgCJiZmYEoivD7/Xm7K5w8eRKRaAyr20GIbApLG2msBQXEnAJcLidCacBfB4BnMdTTiI/vRRBMZpDI/1YPKNUnemqkOGA1nMJGOAWHAzjSksRqOAVOEOGw2cGKNrQ1NSIjppCua0ej14neOje6mjyaARj5OrkTYcxkMrhx44a06LMKURQxPz8Pm82GkydPAihdVKryPPLHtPu69VDhU0GtegvZZF9ZWcHo6KiuS6eQxHQgt/yYfIxi0hQI8qT3xsZGTReLPF2Bczfh2nIckZQT21EXPA4XliJh9LR40VHvg8MeR4ODx3tz2wilOAQZEdGqFb3aghUA0Q6wLHAvyMBm26uRevSQD7ABvCBiPZLCax+vo97tgCCKuDJ8GA8cb9UNwCD1X0nYfiQSwdbWliSMDMOgoaEBu7u7SCQSkih6vd6ylvVaXl5GJpPBwMCAYbEwKoxkr02NeDyedd+W2pVKyYUKnwrK6i08z0tWEikNpkchdT3lgTJqYxQrfKRQNkl639jYUH2P8shNm8OBv/voLkQRSLE8NsIMODaFoZ4WzO5mkAaPw631mA4xYDk7Gnx1SCeSqIVuB7WCFJEriHA7bUhxIjaiafzrzC7q3E54XXYMdPvgctiR4QT8w+0tDB5pRL1be2pR6+QO7C3Mbty4gaGhITQ0NEgWI0nyJq5UEtkoF0QjZcH0uHfvHqLRKM6fP19SocknYAzD4P9n78zj2yjv/P+e0UjyfcSxE8dOYjtOQk47F0c5ChSWhe0vLYWFkHK0EKCUdKHQQti2NIVSlkIPWtiFLYGFciQQCuHMwdVSCAmE+HZ8O75vy7ZuaWZ+fygzSLLkW4Ek+rxefsWOxvOM5Hnm8zzf4/OprKxk6dKlE94xTmW/8ImCKPGFgD/xuVwuysvLycjIwGazjampdbw5PpfLRUlJCampqWEtfyZDfIODgxw6dIiFCxfq4VttXH8EtytY7B5sLpnM5BiKm/sxyC6MsXFY3AIuWUFWVVotLrqtbswG0ffwk6OkdzxCVsHp8f1tRQEkg4DF4UF1qGRbY5meYNJ7/exueUTiCwWv10tRURFz587VfRpH6mXzlwWzWCx0dHTgcDj0xvNQpBhOSLqtrY2enp6AyuyjAY/HQ0lJSVgxChjbjnEkYoxWdIZGlPhCQHNo0ApMtHaJ9vb2Md1I4xWKrq+vZ/78+cPkxyZyzuCmc82uaNmyZQGhp+Dj/NsVtIdDnNmA0SDQ3T+Ix+nAGBMHgkC/3YNRFOkecpOeYMJiF7G6ZPrs3uhe7ziG9re1OGRijR5sbgWb24vL20eMUSQ7JYa50+JIjh1fZacsyxQVFTFnzhyd9EZCOPUT//NppOh0OnW9TKfTqefTNELUJAaXLFlyVD34tPecl5c3YvHOZEOpTz31FDExMVx//fVTcNXHD6LEFwLJycl8/vnnxMXFcf755wcUmIzFPmSsKyxtpTpz5swRSc9/7NHgf1xzc7O+kg3OSfqHWLXjTSZTwIrXZBA5M1Pg5eIhZmdMo7rbxuzUWJ9ZrFFEEkEQBbyyjFdVo6R3gkAFOofcxJsMCAJY7G4EBFweX3TA41UwGXz3UV23jYr2IWJNBk7OSSEpJvA+1AggKytrRK/H8cBgMIyJGDs7O+ns7CQ9PZ3GxkacTqe+Ywy1W5wqTz5FUSgpKSErK2vUeT8aRiLGt99+mxdffJGdO3dOaozjEVHiC4Kqqrz77rtUVVVx1113BYRbtNzdVCTY29vbaW1tZc6cOWNaaY6H+LxeL9XV1Xi93rA5Sa3tQcuZGI3GYQ7U1dXVTDd6+PVlaxhwysiKyoDDQ32Pne0HW+ka8jDo9CIIIgbkqEfeCQRZhUGX7y+uCCqiALEmA4MOL5UdVlbPTaGkdYCnP2nGaBDxyiqfNPRz67l5JBzR/pRlmeLiYjIzM8nMzDxq1671unZ3d3PyyScPU6HRHBa0HWNvb6++g9T6e0O1agS7uIeCqqpUVFSQmpo6TIJxKvHZZ59x3333sWvXrrALgBMZUeLzg6IofPe730VRFL7+9a8Pk/mZCsd0LbRps9koLCzEYrGEdIMIxliJT1VVHA4HaWlpzJ8/f0TRWm1yJyYmBkxYWZYpKysjPj5ed5hIPsL/Q04vr5V0kBpnwiwZaOl3EG82MOSK0t6JClkFWVXpGnQTb3LgkX0h9N0V3STHGok1GrC6vLRanJS1DnBqXpq+68nIyIgoAYTC4OAgVVVVusRgMEZzWPD35NOI0d+s1t/F3Z8UY2JiqK2txWQykXPEMDsSqK+vZ+PGjbzyyiuT3lEer4gSnx9EUeTOO+/E4XDwv//7v8NeHw/xaeovwYRSWVlJTEyMXj02WRshf2hKL5Ik6T6GoaAoCiaTiaSkJFpaWgJCPGazmcHBQVJTU5k2bRpOpxOz2ay/j6pOK/0OD7lpcXgVlTnTYukYcNJj9eCJao+d0LC5ZVp7rXQ0N1Du6qLXYseLQFGvT6Tc4VbYXdnN9AQzW/9ehhQTywXTE8g6igUYVquV8vJyCgsLJyzcPJonX7BZbXd3t25BJMsy8fHxlJWVDSNH/3k2UXR3d3PVVVexZcsWvRcxiuGIEl8QCgsLqaiomLQZrUZo2o3sLz/mv8IdTyHMSNCEpufPn09tbW3Y47S8nslkYsGCBQGvWa1WSkpKmDVrFkajUZ+wGjEajUYq+kX6LC46+630OmQUBSRJZEaiiZaBqI/eiQwVsHoF/nTAQWail3+Zl8C20j6cXhlJAFGForoOiuramRYnEe+RePzvdVx18ixOnZcRcfUSu91OaWnpsEKvqUYoYmxra0NRFAoKCoYV33R3d+NwOHC5XKiqiiRJIfOLoxGjzWZj/fr1/PrXv2bNmjURe3/HA44J4tu5cye33HILsiyzYcMGNm3aFPC6y+Xi6quv5sCBA6SlpbFt27aAUEJTUxOLFy9m8+bN/OQnPxl1vMma0fofazQadYkzrTp0oucMtyrWLJFGUpMB30pUm1jBecqBgQEqKipYunRpSDFgWVHY9mkL/2jvprrPg9srk2CEAZev0TmKKADcskqn1UO3zUOTxY0ggscrkpZoYtGMBCqaezEaJOZkJPusxOxu3vz8MEZLc0BhyUTyZyPB6XRSUlLCkiVLwu7UIoXu7m5aW1tZuXIloijqOfVwLQwejydgxzg0NBQQSjUajQHtGZWVleTm5nL33Xfzve99b8yC/CcyvvLEJ8syN998M3v27CE7O5s1a9awdu1aFi9erB+zZcsWUlNTqa2tZevWrdx5551s27ZNf/3HP/4xF1544ZjH1NoZgjERKbKenh7q6+vDktJkd3zt7e20tbUNE5r2b7vwr9z0b1fQ0N3dTX19PYWFhWFXwp81DfBRg4Xc9ET6HF4ae+z0On1+etEAZxTBUFToc3iJNwp4FJWmXgeNPQ4EwCSBscNJQXYSDoeAxaVySJnO+YvSyUgwDissCZU/G0+Y0O12U1xczMKFC0Mu6iIJi8VCfX09K1euHPOO1mg0hiVGVVV1dRutMnX79u1UVFToBLt161ZycnJ45JFHIrqzPZbxlSe+/fv3k5+fr+es1q1bx44dOwKIb8eOHWzevBmASy+9lI0bN+oP/ldffZW8vLxxVTbFxMSEdDgYj+uCKIq0t7djtVpZsWJFWImz8RKftipWVZWGhgasViuFhYUBk0rLB2q6neHaFcBn8NvR0cHKlStHlGFr7nP4bGhEgTiTEUWNkl4Uo8PmURHwMx4GnF6Fig4rtd1WzJKBk2YmcKDZQknrIHf963zSRsmf+RNjqDChPzEajUbq6urIz88/qrqf4EsdHDp0aFSJw/FAEASdGJOSkkhPT2fx4sXExMSwZcsWBEGgv7+fhoaGqPnsCPjKE19rayuzZ8/Wf87Ozmbfvn1hj5EkieTkZHp7e4mNjeWBBx5gz549PPTQQ+Mee6xmtMFQFIXBwUEkSRpVDWIiTg4Ahw4dwmg0hpRY0ohP8xXTJktwu0JtbS0Oh2NUc1uAGUlm3F4FVVWZHm/Eq0RJL4qxIdx94pZBUWXaBpykxZuwuWXK2gb5+vzpYc8lSRKJiYkjhgk1YrTb7XoVZW1tLbW1tXqYMHjHaDKZprTARjPPXb58eUQJaOvWrXz88ce89tpr+vyeNm1a1Hh2FHzliS9Uj1vwDRrumF/+8pf8+Mc/HndMfzxmtMHwer2UlZVhNBqZOXPmqHmJ8Uw2Tbm+urqa9PR0srOzQx6nOTmoqhqyR09RFMrLyzGbzWPWJjw1N5WK9iHK2gZxeVUkg+/BFcWxgq+onZIKsqxS0T5EblrciFdY1Gxh+8EOnB6Zr81LZe3ymUhB80vbDSUkJFBUVER+fr5eTKZ5S/rvGC0WCw6HQ4/wmEwmnQz9ydFoNI55rmqh1cWLF0e0h+7999/niSeeYPfu3VPqg3gi4CtPfNnZ2TQ3N+s/t7S0DOv70Y7Jzs72JcwHBpg2bRr79u1j+/bt3HHHHVgsFl3Db+PGjaOOK0kSHo8noM9ntKpOTQh67ty5eDyeKanW9IeqqpSVlTFv3jymTw+/KhZFEY/Hg9lsRpKkANLT9AHT09OZM2fOmMd9v6qH4pZB7G6ZwuwkDnUM0mvzfRbRnd+xgK8g6eHrAeyze1CB9AQTy7JC5+Dqum38z4eHSYqRMEkib5d1IYkia5fPHHasoiiUlpaSnp4e8KwQBAGTyYTJZAppbKuqKm63RVtrVAAAIABJREFUO2DHqOUYNUk/s9kcQIzBBrxer5fi4mLmz58/YfPcsaC0tJSf/exnvP3222F3v1GEx1ee+NasWUNNTQ0NDQ1kZWWxdetWnn/++YBj1q5dy9NPP81pp53G9u3bOffccxEEgQ8//FA/ZvPmzSQkJIyJ9MDnyTc4OBhAMCMVt1gsFqqrqznppJNISkqio6NjTE7oY8XAwAADAwPMnz9/RNKTZZnk5GSqqqr0ie6vXN/c3ExOTs64moaLWgZ48UAbGUkm0uKNfN48QKJZYsDhDajoNIBuPxslwyjGAlEAsyQSYxQxGEQGnV5S44Y3lR/qsAIQf0T1ZVq8ic8OW4YRn6aMkpSUFJAiGQsEQcBsNmM2m0PqZ6qqOsxnUMsxagtil8tFYmIiNpsNRVEiYqfU3NzM9ddfz7Zt246q4s3xhK888UmSxCOPPMIFF1yALMtce+21LFmyhLvvvpvVq1ezdu1arrvuOq666iq9XWDr1q2THldzaAgmvlA7vo6ODlpaWgLi+ZIkYbfbxzzeSOLXXV1dHD58mPT09BHzBdqqNDs7m5ycnIAVbF9fH4cPHyYpKYm2tjaampoQBCEgz+H/5b9LrOq0YpIEjEf0F1PjjFR1WokzS9hdX5BfvFlg0BWlvCjGDlkF6xHVn/puG6+XdPCDs3IAGHB4iTcZiDEaiDP7/P40uLwK0+IDCVJVVaqqqjCbzeTm5k75tWrzJSYmZlihjKqq+i4zNTV1mM+glmsPzi2O5BoRCv39/axfv55HH32URYsWTfl7PFHwlSc+gIsuuoiLLroo4P/uuece/fuYmBheeumlEc+hVX2OFSkpKcOa2IMLUYIrK/1XdRMpWgleFWrGtH19faxYsYKmpqaQ5wzXrqCtYK1WK11dXaxZsyZAe1QzBNVWsP4q9v6ahO4hhSG7kySTT2nG5ZWJMYoICBhijAw4fKGqKOlFMRGowJBLZsgl89ynrdjdCl1DTvodXgyCwPdPm83JOSn8vaaXNosDEDBLIpeuCNzt1NXVoaoq+fn5R/f6jxBuXFzciGopwXZKwfMtlDi2fw+j0+nkyiuv5K677uLMM888iu/w+MMxQXxfBrRQpz/8d2Sa/Fi4ApGJmNH6E5+iKNTU1KAoCsuXL0cUxZBkOlq7Qltbm948G6xLGM4QVINWOp6YaqOsq4WmniFURcEkyEw3CcQZRZyKgNcr4PCqCOoXoU5JjDa2RzF+DDi8bD/QwqKZ8WSmxOFR4ImPm7hv7UnccX4+Ja2DuL0KC2bEMzPpi+iH5q6wZMmSo+4/19jYiKIoo0qEjcdOyeFw6OLYr776Ktu3b0eWZaZPn05JSQlDQ0MsW7aM1atXR+ItHfeIEl8YBLuwa9DCh6WlpcycOZOsrKyQvz8ReTMNXq+X8vJykpOTmTt3rj6Rg4/TSC9cu0J9fT1DQ0Pjap71h3/p+K8vyeBQpxVZUZmXHse+hn6e299CgqoQ53Dg8HoxiuBUQMTPuTuKKMYJu0dFUj309fWhqioDTpUPD5RRkJVEVmwsMQkxxIpe3VuvubkZi8XC8uXLjzrptbS0TJlzezg7pWXLluFyuVAUheuuu47Dhw/T0NBASUlJlPgmiCjxhUE49RbNPyyU/Jg/JqLyAr7keGlpKdnZ2cycGZi49yfTkUhPURQqKysxGAwUFBRMycPAJIks96u4O++kdGYmxVDdZSXeZODpT5rxKgqHex1Y3T7W829cjiKKsUIQIC4+gcQYI15ZwWt1s3zBXNJiVD13pjmuu1wuvF4vaWlp1NfXT6nM2Wjo6uqio6ODFStWRHScxx57jPb2dp577jkMBgPLly+P2FgnCqLEFwYpKSl0dnYG/F9vby8ul4uVK1eOKn00kVCnpum5YMGCkFVlBoNBX/lpvoAGgyFg0nm9XkpKSpg2bVrAbnGqIQgCy7KS9PJzoyTyk5fLcXt9VCcKPtmqKKIYLyQRPm8aQFEh1iRy2cosFmYNX2RqRV8FBQV6EddoMmf+X2azecLzo6+vj8bGxjGJP0wGr776Km+++SZvvfXWlI5z7bXX8sYbb5CRkUFZWdmUnfdYQZT4wiApKYmamhr95+bmZrq7u0lOTg7QxAwHTVZsLDAYDFgsFrq7u0cUmtZ2fFo+MLhp1eVyUVxczJw5c4btFiONoqZ+BFUl0STgkH2h1ijxRTEReBVf1SaAWRH5Z10f3zgpnfkZX8yL3t5enXiMRiMmkymsUEWw6PPg4KC+WwQC1Fz8v8I1rQ8NDVFdXT2iFOFU4OOPP+YPf/gDu3fvnnL1l+9973ts3LiRq6++ekrPe6wgSnxhoDk0aEUmsixTWFhIZWUlXq93TOQ3VthsNmw2W8gCFH8IgkBPT4+eJPdfudpsNsrKyli4cOFR1yR0uVzUNnViEEXiY4x4HR7snijrRTF+CPiiBUaDiMur0GV10211c+X/HeDcBdPZeHYe8TipqakZVV9Ww2iiz/5qLlrrj9Pp1Ptw/ZvWRVGkqamJpUuXTukzIBiHDh3itttu44033ojIfD7rrLNobGyc8vMeK4gSXxgkJydjsVh45513WLBggR42nAoXdg1aAYrD4SArKyss6Wn5vMTERJYuXao30WorV5vNhsvlIikpic7OTgYHBwNWrlPZPBsMu91OSUkJ/7I8i4Pvt2JzefXVehRRjBcqvmrgAac34P+GnDLFLYPc80YFV+a6WLNqhT5fmvsdvr5So4HVc1MwSWPPt41FzUWbb0NDQ9TX1+vRIK/XO6wX1v/7iYYmOzo6+P73v88zzzwzZnWlKMaHKPGFgdVq5Z///CfnnntugLffeKo1IXxjur8b++zZs/Xm82AEF7EEhzw6OjpoamrSq7v8VSW6urr05ll/n7PgHqGJ5jkGBwcpLy9n6dKlnJKQQL9H4qmPmwAZSYxWdh6/OPq6nwIw4PBgVJxknPOFUERxywAP7K7FeySuPj89nl/+28Jxkd+I4x4hNkmSqKmpYdmyZQFFbcG9sP6FNxPxFxwaGmL9+vU89NBDFBYWTsl7iGI4osQXAnv37mXjxo2kp6dz2WWXBbw2nmpNf3sgf7jdbsrKypgxYwZZWVm67FEwRmtXOHz4MP39/axcuVLf1YVSlQAf0fqHc7QxNZslLc8RFxc3LM8RCr29vdTU1AR4+F37tTnkp8fzm51VWBwydrc3Sn7HJY6+7qcK9Du8eGSRLfvauGKNyOLMRLZ83HQkNOq7pqpOK/sb+zkjP23KxpZlmeLiYubOnTusknu0Xlj/eTdS4c0LL7zAjBkzeOONN7j88ss599xzp+z6oxiOKPGFQHl5OTt27BhGejAxF3Z/4rPb7ZSXl5OXl0daWpp+XKjGdK/XG9JdQVVVDh06hKqqFBQUjKmU2mAwkBDG58w/z2G327HZbPT09OgahNqqVyNFrbm2sLBw2A701NxUFmUm8c/aPiSDiFeJMl8UUwerW+GNsm52VnQzb3o8Hlmlz+5GVgBUREGgzx46ejIRaMLwGRkZzJgxY9y/P9K8gy9slBYtWsQrr7wCwHvvvcezzz6Loig8/vjjnHzyyZN6D1EMR5T4QmDDhg2oqhpyZ6dZA40FWlhUS4JrQtaLFi0KSLQHN6bLsqxXbga3K8iyTGlpKUlJSeTm5k5Ju8JoeQ7/cI7WsJuQkEBxcTGKogyrirvrnCx2ZyfyfnUfblnBI6vEmw20DbhwuGW6rFMn3h3FiQcBXxi9vseOrKqIAiTFGJEVFZes0Dk4NrPo0aAtMBMSEsYteD1WGI1GJEmisbGRnJwcHn30UX2+q6o65srw8eKKK67ggw8+oKenh+zsbH71q19x3XXXRWSsryKOGeLbuXMnt9xyC7Iss2HDBjZt2hTwusvl4uqrr+bAgQOkpaWxbds2cnJy2LNnD5s2bcLtdmMymXjwwQfHFUYIZUY7VvFp/51cZ2cnzc3NIY0p/YlPlmUURQnZrqD5fGVlZY3LXWGy0PIUra2tiKLI6aefPszqyD+Mau/vZ5HRQf4C38pby2181mnmuWILMZKAy6tGm9ujmBC0+0YUVN9OTwCPLGM0GJgVL6FOUSi2vr4egLy8vCk5Xzg8/fTTlJeX8/LLLwfMK0EQJrWwbW5u5uqrr6ajowNRFLnhhhu45ZZbAHjhhRcmfd3HMo4J4pNlmZtvvpk9e/aQnZ3NmjVrWLt2LYsXL9aP2bJlC6mpqdTW1rJ161buvPNOtm3bxvTp03n99deZNWsWZWVlXHDBBbS2to46pnbDBRPfeIpbNE+/xsZGBgYGhglZ+x/n9Xr1Ahd/oWkNNpuN0tJS5s+fr4dIjxY0JRhJkli6dOmwyaiVi4dq6lcURa+KOzvBQVY8/PafPQw4ZPqjG78oJghf24OAWfKR39xEAVXxMuT2YrK2U1FhG3NfXig0Nzdjs9mmRIpsJOzatYsXXniBXbt2TXn1tSRJ/O53v2PlypUMDQ2xatUqzj///IDn5omKY4L49u/fT35+vr7yWrduHTt27Aj4A+7YsUN3YLj00kvZuHEjqqqyYsUK/ZglS5bgdDpxuVxj6sGJiYnB6XQGJK7Hk+PTen5iY2NZtmxZ2FycIAi6Y7rJZBpGehaLhcrKSpYuXXrUTSdlWaakpITU1NQJKcH4V7UBZGVlkZlt5eevH2Kow6pX40URxWjwVwNSAZdXZXqCicQYCdkg4EbFaIY3mwRaZZX1y816X151p5W2IQ8pZoGT0mP0Plj/SkuNeDo6Oujq6mLFihURJb0DBw5wzz33sGvXrrDFMZNBZmam7teXmJjIokWLaG1tjRIfxwjxtba2BsTYs7Oz2bdvX9hjJEkiOTmZ3t7eAD+9l19+mRUrVoy58VQTqg4mvrFUdXq9Xrq6ukhISGDhwoVhJ5BypPjDbDZTXFysf69NRrfbTXd3NwUFBWEVXSKFSIVWF85I4IrVs/jt7jqGnN5oyDOKsUEFoyiQYBZxyyqyomL3yNxx/jwWZyZy28vlSKLPsuiTJitGo4k7/mU+rxV38H9VVlSMoKr8qzmJy3Om4XQ6h7X9KIqC2+0mMzOTtra2iOl+NjQ0cPPNN/Pyyy+TkZExZecNh8bGRg4ePMgpp5wS8bGOBRwTxBcqwRtMJKMdU15ezp133snu3bvHPK5GfP4ux2PZ8TmdTsrKykhJSSEhIWFE0tMqNzXrIX/z2ObmZgYGBkhKSqKiomJYP55/68FkdAdDweFwUFxczLx580hPT5+y82ros3pwemQMAniD/nRmScDrVRlb00gUJwKMom8+y4qK06vilVUQINlkoL7XjlEyIKuQfMShPTXWyCcNFmwuL0/ubSLBbMBoEJEVlZ2H+vhmQRazZwW2JgwMDFBRUUFBQYHehhDcfiBJUkjPvPHMv56eHq666ir+8pe/MH/+/Cn/rIJhtVq55JJL+OMf/ziqxvCJgmOC+LKzs2lubtZ/bmlpGbYD0Y7Jzs7G6/UyMDCg99y0tLRw8cUX88wzz4zqmeWP5OTkYQ4No+X4hoaGqKysZOHChbhcrrCFMBrpGQwGJEnSV5NahWVjYyOCIPC1r31tWFWnfyGJ5vCs6Q767xb9yXE8KhJWq5XS0lIWLVoUUix7KuBVfGGqIacXq0v2K1iAGAOYjdDjAFEUouHQKPAoPvJTVHB4FN35o2vIjdOrEGcyoKpf5OQ9skKC2cCQy4tXUemz+fLnKbESBlFgyBk4h202GxUVFQF9qSGvYxTdT5PJFFb3E3ztTN/97nfZvHnzUdl9eTweLrnkEr773e/yne98J+LjHSs4JohvzZo11NTU0NDQQFZWFlu3buX5558POGbt2rU8/fTTnHbaaWzfvp1zzz0XQRCwWCz827/9G/fffz+nn376uMYN5ck3kvh0T08PDQ0NLF26lLi4OHp7e0OGRUdrVygvLycuLo4FCxYMW0WO1o+n7RbtdntY9ZbgJnX/1arFYuHQoUMsX748oqHVBTPiiTUZyEqOQVUVmvqdzEmNJV5w0mVXSU6I58IZCbRYHHx+2EK/82jt/46+KkkUY4PHbwEkCEdsr1Tot3lYMzeFhTPiqeq0ofXz/ejsXJwemc5BF26vgihC+6BAzrRYZqd+QW5Op5PS0lKWLl06IunB6Lqf2vzTpAS1fthnnnmGDz/8EIA5c+bQ1NTEm2++yaJFiyJWNaqqKtdddx2LFi3itttui8gYxyqOCeKTJIlHHnmECy64AFmWufbaa1myZAl33303q1evZu3atVx33XVcddVVuk/e1q1bAXjkkUeora3l3nvv5d577wVg9+7dY4qrhzOjDYWWlha6urooLCzUV3ehwqIjtSt4PB6Ki4uZOXMm2dnZYxrXH4IgYDabMZvNIXdqo+0WwdcWkpWVhc1mQ1XVSWkOjoRzFkynqtPKzvJuEOCU3FQunu1hRtoscnJyAgh/Z3kHt71cOeXXEBpR0jsWIIkCogCyAmajiNEgsH5NFnvrLSTFSpySk8r8jHj+8G4d0+KMDLq8uDwKggDzM+JJjPE9+rQ5t3DhwkkXjo00/1auXMntt9+O2WzmW9/6Fo2Njezdu5eGhgY2btw4qXHD4aOPPuKvf/0ry5Yt0+XPfvOb33DRRRdFZLxjCcIoDZIndIzpz3/+M7IsD2vs3L9/v66moKoqdXV1uFwuTjrppACSsFqtHD58mCVLlgCM2K4Q6ZzaaGhpaaGtrY3c3NxhfXnabjF4pzgVucVBpweHy0NLbSUZGRkhG4VlReXrv/8nffapEQeP4viAySAQbzLw2Prl7Czv5v3qHkRRQBBg878tZNWcFO57u5p9jf0kxRhRUbG5ZBZnJvKbby1ClmUOHjzInDlzIlpgoqoqf/zjH6muruapp56asiIZp9PJWWedpZvxXnrppfzqV7+aknMfJwj7YDomdnxfFlJSUkJad2ganKqqUlFRQVxcHIsXLx5GAFoFqJbP0/J3wTe+llRfvHhxSOWUSEJVVRoaGhgcHGTVqlVhd3fjyS36E+Rou8VYA1RVl5OdnR1QROQPj6xQmJ3Me9W9k3uzURw3EICTZsRzy7l5yAq8X91DQoyEKAg4PDL/tauWFzes4usL0viovg+HR0YA3LLCuQunoygKpaWlZGZmRryq8sUXX+TDDz/ktddem9LKULPZzHvvvUdCQgIej4czzjiDCy+8kFNPPXXKxjheESW+ERAu1Kmpt1RVVZGZmRm21F8LdYbT3ATo7u6mrq5u1KR6JKCqKlVVVSiKoleVhsNouUWtSV1L+GvEONJuEaC4uJi8vLwRd7nP7GuhuGVAL2iI4sSGAKTGwLo8mXhrKyW9vvYhRRZAFIiRBPrtHmRV5Yx5adx6jsy2A22oKlx96mzOXZBGZWUlycnJZGVlRfRa//73v/P444+zZ8+eEb02JwJBEPT56PF48Hg8Ee07PJ4QJb4RkJKSMqyqE74Qrl2wYMEwtXZ/CIKA0+nEYrGQmJg4LKfX0tJCR0fHqAa0kYCiKJSVlREfH09eXt6kJowmYj2aM4TdbteJ0Wq1MjAwgNlspqWlhd7e3mHkqO0Wy1oHGHDKI5KeCETlsE8MxBhF8jISueCM5SgeF3KcBYoHcbjciChY3Sqz4gWKPv+c2NhY8mNj+fV5M/T7qq6uDkmSAuzGIoGysjI2bdrE22+/HTHhCVmWWbVqFbW1tdx8883RPr0xIkp8IyBUO0N/fz9DQ0MsXLhwRNLTLIXy8vIYGhrSbYC0what6GXu3Lk4HA7deuhowOv1UlxcHDanNtUI3i1q7RJr1qwhMTEx5G7Rbrfrlk5mrxflSC46eNcnHWl3mJVionvIg0uO7gmPV2h/+7R4IxvPziXWJIFJYkV+PP95kZnfv1uH26uycFYsv/rmSUyLFQPC8xaLhb6+PtxuNzExMRQVFYVsPZgK6bDW1lY2bNjA1q1bI6qrazAYKCoqwmKxcPHFF1NWVsbSpUsjNt7xgijxjYCUlJSAUGdHRwctLS1kZGSMuEPzb1fIzMwMCCFqO62YmBimT5+O0+mkv78/wAIoVB/eZAxj/eFyuXRvsYnYrEwWWj7Tv11itN1iZs4gf99SgtWtBJCeUYD0OOixQ4slKvx5oqDV4uJH20rZ/M2FXLjEdw+fvWA6Z81Pw+GRiTd98Vjz15Btb2/H4XBw2mmnIQhCwIIrVOuPZsXlb8llNptHzdNZLBauuOIK/vznPx81ebCUlBTOPvtsdu7cGSW+MSBKfCMgNTWVwcFBVFWlsbGRoaEhCgsLaW5uDtvEPlq7QklJCenp6cyZMyfk7yuKMqJhrH+DrH9ocCyrVLvdTklJyagh2kghlHntaDAYDGSnp/LyjWu49aVyarpsKKqKQRDISDTj8sqoTJ3/WhRfXWiLHskgYHfL3PtWNd84KR2TwUdEoiAEkJ4/enp6aGlpYcWKFTpxjbbgcjqd2O12nE4n3d3dOJ1OXcEl2IpLy+OnpqZy5ZVXcscdd/D1r389Ip+Dhu7uboxGIykpKTgcDt555x3uvPPOiI55vCBKfCMgPj4eh8PBq6++ypIlS3Sl9nCyZSO1KzidToqLi8nJyRlxpyWKIvHx8SGbxyfToD40NER5efmXInQN0NXVRWNj44TzmbNT43j5hjV4j1TTdgy66RhwsmVvE3vr+pFEn0q/9nA0GwTcsooIUemz4wjCkS9RFHB4FGwuL6a4ke+ngYEBamtrWbly5ZjDmAaDYcR5GNzy88knn/DYY4/R2dmJqqo8//zz7N27l7y8PDZu3BiRopP29nauueYafbF92WWX8c1vfnPKxzkeESW+ETA4OEhvby/V1dV8+9vf1m/eYNmy0doVhoaGKCsrm7QE2Fgb1P2LSOx2O3a7HZfLRXJycoDwrhbGiUSDuj/a2tpoa2tjxYoVk85jSkc+29mpPvWNV4o7EAQBVQVJVPEovgdjTlocPTY3yTESh/scRFN/xw9UVUVWISlOIjl25PvJarXqUmRTVUAWyrg5JyeHiooKFEXhwQcfpLW1lYaGBlpaWiJWabl8+XIOHjwYkXMf74g2sIfB4cOH+fd//3e6u7spLi4OuHm7urqw2Wzk5uYGCE2HalfQwnvLli076u4K4MtLNjU1UVBQgKqqOilqBOl0OvXQbKjc4ng8zELh8OHD9PX1sXz58ogQbEX7EDdvK6VnyOUTGxNg2awktlxVyA3PFdPYPYTFqUSJ70vD1EnAiX6nMUsij68vYPXc8AtJp9PJwYMHWbZsWcg2nKnE448/zt69e3nhhRcicp/Lsszq1avJysrijTfemPLzH6eINrCPFzU1NTzyyCPceOONenWhhuDG9GChaQ1tbW20trZ+Ke0KAE1NTfT09ASEeILd3zX4h27sdrtecKP1Bmlq9P5h1JGsWlRVpb6+HpvNRkFBwZQ27vpjcWYif/luAa8UtdMx6OLU3FTWLp+BySBy3VIjD+8zYPf6hI2DYZYEvLLPCV4U4D/OzqW6y8aHNT04ZBV3sG1EFBPA5EjPIIDR4JMkM0sGjJLA9afP4V8XzyAlLvxuT7PUWrRoUcRJ77XXXmPHjh3s3LkzYtGThx9+mEWLFoVsr4pi/IgSXxicd955ACQkJGC1WgMUVTRn9XCk56+GsnLlyoiHEoOhqiq1tbU4nU4KCwvHRDqjuahrivR2uz3AqkVRFL3gxr8CtbW1FUEQIu5gDZCfHs9Pz8/Xf9YUdbKSY3j5piW0WJzc8UoF9T12PLJCSqyR6QkmOodcxEiiT93faOCSlbP4/bv1SJJEeqxIn80dkjCjODpIMkGKWcDiUokzCgiCitOtUNpi4bKV4VsEZFnWhREi5S6i4ZNPPuGhhx5i9+7dYReVk0VLSwtvvvkmP/vZz/j9738fkTFONBxTxLdz505uueUWZFlmw4YNbNq0KeB1l8vF1VdfzYEDB0hLS2Pbtm16k+r999/Pli1bMBgM/OlPf+KCCy4Y05iaeos/8YmiyMDAgG40GxcXp5c5K4rCoUOHEASBgoKCo66koCgKlZWVSJLE0qVLp2R8TXklLi6OtLS0gNf8E/12ux2r1UpNTQ2KomA0Gvn0008DysE1gpxq/0ANoRrz50yL4+lrVlDZbkUQYNHMRGxuL/fvqqGsbYi8lFj+84L5pMaZWDM3hX/U9IIKyTESDk+0TeLLQv6MZNITjXzaaDniiqKgKDLvHOrh0/p/kJ8isn5JHNMSA30p6+vrmTVrVsQ1b6uqqrj11lt5/fXXI1olfeutt/Lb3/6WoaGhiI1xouGYIT5Zlrn55pvZs2cP2dnZrFmzhrVr1wb0yWzZsoXU1FRqa2vZunUrd955J9u2baOiooKtW7dSXl5OW1sb5513HtXV1WPaiSUlJQX08smyjNlsZsGCBTidTnp7e2lubtbLnF0uF3FxcWRkZNDd3a0/6I/Grk+WZUpKSkhNTWXu3LlHhXT9E/0JCQmUlpaSnZ2tLziCNT67u7v1YhuYuMZnKGj6i8nJycNUOcySgcLZXyxeTJKJ3168ZNg5/t+yGTT3O9h2oA23V9F3hFFMLUaTn4uVBB78zmI+rO3lQNMggiigKgaG3F5MBhGnaqC4V0VoMHL3eZk0dA8iWgc53NaJUQS3u4mWlhb9/goO0U92PnZ0dPD973+fp59+mrlz507qXCPhjTfeICMjg1WrVvHBBx9EbJwTDccM8e3fv5/8/Hzdu2rdunXs2LEjgPh27NjB5s2bAbj00kvZuHEjqqqyY8cO1q1bh9lsJjc3l/z8fPbv389pp5026rj+ep3+7QrBq0mtMXzevHkkJiYGtBv4q5AEP+Tj4uImXUACX+Q0srKyIqoUEQ6aGsyMGTMCLJXGo/E5MDBAR0eHrnAT6vOKjY3FZDIN+7w00p8+ffqk1GgEQeDmr+fygzNz6LG5+ff//QyP7A4ojhHwFdGgRmXSJoI4k4FYo4hXURlwhHHcEAS6hlxcsmIWjb0QzrxyAAAgAElEQVQOXivtwONVMBlEkmMlBEHAIKp8dniAK5+vxOlW6Hd4MBkEYk0SZ+Ync883FyJ73AF567a2Nj1EP577yx9Wq5X169fzwAMPsGLFigh9Sj589NFHvPbaa7z11ls4nU4GBwe58sorefbZZyM67vGOY4b4WltbAx5o2dnZ7Nu3L+wxkiSRnJxMb28vra2tAYrl2dnZtLa2jmlcTa9TI71Q7QpWq3WYdmeo3ILX6w1oN7BYLNjtdv3coSbhSAUkGr5sSyO3201RUdG41WBG0/jUPi/tM7NYLDgcDtxuX/hRK7iJiYmhra2NGTNmTJnosEEUkBUVr6IQbABvEAWMBgGvrKBEmW/ccLhlnG4Zo0FAAJJiJNyygsOjYBAgKdaIoqj87z8P8+i65dzxL/nccm4unzRYuPv1Q/p5VMDm9rkuOD0yqgpOr0qMUeXD2j7eKOvkO4WziI2NDRmK9L+/tIWX3W7X7y//aITH46Gvr4+8vDxuvPFGbrrppjGnSyaD+++/n/vvvx+ADz74gIceeihKelOAY4b4QrVdBK/Kwh0zlt8Nh6SkJPr7+3UtzWAS6u/vp6qqiqVLl45aPSZJEomJiSEbyMMVkDgcDiC0YktcXBwOh4PS0tJJ9whOFFpj/rx585g+ffqUnnu0z8vlcjE0NERNTQ1xcXFYrVY+/fRTPb8YnFfU2jPGip3lXXgVlWnxRhTFZ2mTEivxtx+s4Tc7a3iztItQATuTQcCjqMxIMNIz5CHqIvgFBMAkibi9Cl5FxSQJmCURl+xbQUgGEZNBxK0qOP0Ki8ySgdPyUsmdHkdNlxVZ8em0+r7Qd+SCAIoKqqoccWMPj5Hur+BoRE1NDY899hg1NTW6wPquXbvIy8vjpz/96VG3E4ticjhmiC87O5vm5mb955aWlmEhPe2Y7OxsvF4vAwMDTJs2bUy/Gw5dXV38z//8D8888wy5ubnk5eXp/+7fv5+WlhY2b9486Yqu0QpIQim2DA4O4nA4iI+PH+ZwMFUh1JFgs9koLS3lpJNOOuqkK4oiBoOBxsZG5s+fP8xTzePxBPQsaguJ8eiheo9s5wRBwGAAIwJGSSTRbOSO8/M5cHiAVotT71RTgVijgNFgIMUo8tQ1K7C5ZG7bXkabxYnWHXG82ytNjzOSEGOgsc857DVRgASzAY8kkJsWR5/di9XpwSD4HNXNBhGP7HNKX1sQGD0wGUT+54rl/O1gO+0DTgpnJ/O73dX0Wt1Iok+pRxtDQCA/feJ9s8HRiMzMTD799FOys7N57LHHcLlcNDY2Ul9fH7FqzmCcffbZnH322UdlrOMdx0wDu9frZcGCBbz77rtkZWWxZs0ann/+ed3dHODRRx+ltLSUxx57jK1bt/K3v/2NF198kfLyctavX8/+/ftpa2vjG9/4BjU1NWNOcKuqysDAAHV1ddTW1lJbW8trr71Ge3s706dPRxRF5syZQ25urk6KeXl5ZGVlIYpiRMinq6uLhoYGCgoKkCQpICSo/RscQg1+yE+mt05To/myJNBcLhdFRUUT2mkG66Fqn1mwHmqvR2LT7g7cXhVRBFEQueGMOXzvtDkoisKBzw/yepPAR01OJIPA0sxEBpxeEmMkfnBmDgtnJFDdaWVneReqAJcUzuQP7zXwYU0vblnRH9THCrKSTSSZRGp6viDxUA3qi2fE8Z8X5HPD1nLs7tCCcQkmkQ2n53Dpykzer+7Bq6gIqLxS3ImiqFy2ahbfKcwcce709/fzzmeHeLxCxeVRsDg8mI0iMUYDp+Sk8sDFi3Sln8nimWee4c033+SVV16ZEvcGf+Tk5JCYmKi3Rn322WdTev4TGGFvnmOG+ADeeustbr31VmRZ5tprr+VnP/sZd999N6tXr2bt2rU4nU6uuuoqDh48yLRp09i6dateDHPffffx5JNPIkkSf/zjH7nwwgsnfB2///3vKS8v57HHHkOSJFwuFw0NDTox1tXVUV9fT1tbG4qiMHPmTPLy8sjJyWHevHk6QU60pF/z8SsoKBg1dBccQvV/4KuqislkCtgljkX02mKxcOjQoS9NjcbpdFJUVBQRse3g3XVZ6wDPF/dhdcmszhA4Z46PFIeGhnQjU62MPvhvWdQ8wE1bS3B6FBRVJUYS8SoKpiMLLpdXweVVwk4yg8BRV5wJ5WsoiXD+Sencf/Fi6rpt3PR8Cb02d8hrMxsEblqdxOetdva1e3z5txDcF2MU+fC204k3T4xENO3ZwsJCZMEnSxdrNODyKpglkZy02ClbcO7Zs4cHH3yQXbt2ReR+z8nJ4bPPPpvyVEEUxwnxfVUwODhIYmLimCaWLMu0trYGkGJdXR2HDx/W9TP9d4ravykpKcPOrygKdXV12Gw2li1bNumSbK0Hz3+XqJGi1pwfnCdzOBy6BNrRCvH4Qyvk+TLCq/BF9azWv6l9bv4FEdpn9vN3Oihtt+M+og6jIcFkQDKIPpKVfeQXXEBzam4K/71uOac/9M9xNdGPRJYmgy/IajSImCQDsqLg9MjMSIyhc8iJWTIw5PqCpbS90vmL0vnDv39hdVPV1M6rn9azo17xXb/HR9550+P40dm5nL8og02vVLDnUA+SQWDIOTzLGSfBz08xk5UUmIsdS5je4XBQVFQUYG0VKRw8eJCNGzeyc+fOiNl4RYkvYogS31cRqqrS19enh081Uqyvr2dwcBCTycTcuXPJzc0lJydH99q66667pjzcEgqyLAeQYm9vL4ODg3qzvn9j+lSFUEeCzWajpKSExYsXfynFBJoiSEZGRkDLhgatIEL7zDbuaKS6xx0wiQTAKPrcBSRRZFaymT6Hlx7rF43yJoPAzKQYYowiKbESB1sG8Yxh62cQBeJNBoac3pAT1yz5/i4ur0Ki2YBBFLjhjLl877Q5NPU5uPWlUg4FFYQYDQJXn5LNrefO4/3qHuo7+klw9XDpOavpcyq8c6gHVYVzF6aRlfKF1dSBJgsbt5XikRUcQT6KouBradj9o9OIN6Lb//hHI4Ird7X7S5IkDh06xKJFiyJ+DzQ2NrJu3Tq2b9/OggULIjZObm4uqampCILAjTfeyA033BCxsU4wRInvWIOqqjgcDurr6zl06BD33nuv3g/X2dmJIAjMmjVr2G5x7ty5o/YhTQQtLS10dnbqOUWtqjJ4t+jvVxaqqnKihK25tn9ZOUWN9MbTMvHc/hbu31UTsJuTRFiTnUBeikSsQeaUDHC63Pyl1E2VRdWLNFTVR2QGg0DOtDhqu21hc4IigOAjNlVFd6v3+O00BQFiJF++2SMrpMYa+d/vFrBgRgKvl3Rwz1vV2ELk42KMIi9uWMXD7zWwr6EPt6xgNBj48TfyWL9mOPn745OGfp7a24Tbq9A24KR9wIV4RHvz2q/N5qazckf8fVVV9TC9w+HAZrPR1tam234FRySmsqirt7eXb3/72zz66KMBrVCRQFtbG7NmzaKrq4vzzz+fP//5z5x11lkRHfMEQZT4jmX88pe/ZObMmdx0003AEVsWWaapqWlYXrGpqQm3201aWlrIEOpYQ7T+aGxsxGKxjDm8OlIIVfMNDEWK4Qh7cHCQ8vLyoxLaCoWJkB74Pocrn/qcohafsLBkgBjJwD3/7yQuWBxYharlYp/6uIm/7OtAFHy/r6pgNsD9X5N4+KCX2sHASRljFLn3mws5aWYCr5d20mvzcHpeKvsaLeyu7EYA8tPjKG4dwnDE3kBVVcySgY9/ega13Tau2HIAj6zoRClwpDJSgBvOyOGUnBRu3laKyy1jkHwhWlEQ+eSOM1CBh/bUsruym3iTxE/Pn8fZC0KH7P66r5k/vFuPeOQ6fnHhfL5VkDmmz1JRFIqKisjMzCQz0/c7Ho9nXLtF7WssfbEXX3wxt99+O9/61rfGdH1Thc2bN5OQkMBPfvKTozrucYoo8R3LUBRlXOFDRVHo6emhpqYmIHyq5QdjY2PJyckZVoWakZERMI6iKBw8eBCz2czixYunLIQZ7Buo/ev/wNIeVoqi0NrayvLlyyOush/uWrUH7kQUcZwemZ+/doj3q3sQBYENp8/hhjPCy8k9u7+FP7xbhzZnFVUlwSTw8LkJWEzTueP1Btxen+mgqoBBhAe/HseM5Nhhiwmt4OZAk4UfvlCCR1H1J8EZ89J4+LKlvFbSwX07q/F4VRTV16sIkGg2cNJMn/PFrrI27nmrBhkR4QghIwh88OOv8af3G3i1uB2P10eaZklky1WFLM8KFDt3exW+/vuPsLtlRFFAUVVMBpFXf7AmIEQaCqqqUlpaSkpKCnPmzBnT565FJIIXXuFsuHyLDJWMjAy+973vcf755/PDH/5wTGNNBjabDUVRSExMxGazcf7553P33Xfzr//6rxEf+wRA1JboWMZ4CUcURTIyMsjIyOD0008PeE1VVWw2m06ItbW1vPTSS9TV1dHd3Y3BYCA7O5u5c+dSXl5OVlYWP/3pT5FlGUEQpiSEOpKMmX8Itaenh46ODhITEykvLw/pBDHZEOpImCzpAcQYDTx0yRK8ioJhDJ/fhUsyeOKjwww4vCiqilEUWDvPRGFhIZ81DSAZRD10qgoqouhzwEiQvmjRCKWHes3yBP6veAiPrLJydiK/+rf5AMxMMqOqHNnF+RRpJFHg95cu4ZTcVLxuNwZLE4LBgOJVEAEVgbmpsSSaJXZXduPxqr57A1/+8IPqnmHE1211IauqvtsTBd84DT32EYlPVVWqqqqIi4sbM+mBbw5opBaq8jfYQf3AgQP8+c9/pqurS//c6urqyMvL45prronYoquzs5OLL74Y8LVsrV+/Pkp6RwFR4jvBIAgCCQkJFBQUUFBQEPCaqqp4vV5qamr44Q9/qDfT/+IXv6C5uRlZlklPT9dDp/5h1Li4uClzgoiNjcVms2GxWDj11FMxm8369fmHUG02G93d3QEh1FCkOJGcp6Y9OmvWLD20NhmMtZ8sLd7EtutW88y+Zlp7BliS7OHaC9YgiiLLs5KYFmekY9CFrKgYDQIF2cnMSI7V/67B0ApucnIcfGu5HZvdjsvppLq8WI8krJ4hsb/djSiAAYFff3M+X8ubhsfjobi4mK8VLua/c2HTqxX02TwsmpnA7y9d4hMCMBqw4NGX1lqBTaj3Bb4drHhETcmrqGSnjrzba2hoQFEU5s2bN6bPb6wItuGaO3cura2tlJeXs2XLFtra2qivr6e+vj5ixVoAeXl5FBcXR+z8UYRGNNQZxTC88sorHD58mFtvvTXg/xVFobOzMyCEWldXR0NDg64go5FhTk6OHkJNS0sb18Ojq6uLxsZGCgsLx2XgGxxCDW41CM75xMXFhaxC9Xq9FBUVkZWVNSWkNxG0t7fT1tZGYWFhQF61z+bmoXfqaOixU5CdxH+ck0dcCKIZK7xeL3a7nc8a+2jrt5EVp5Bs8Oi77sTERJKTk8Pmyd451M2mVytxeRUkUfBJut24htS44X+396p8xxoEn6Tbj87O4ZpTw+/iNDWiZcuWRZR8AF566SWeffZZ3nzzzSk3jbZYLGzYsIGysjIEQeDJJ58ck0B+FJNGNMcXRWShqiqDg4MBrRlaXrGvrw+j0cjs2bMDCm3y8vKYNWsWBoNB35F98sknGAwGCgsLx6WrOZbrCy6EsNvtw8x0Y2Ji6OjoIDMzk9mzZx+VtpFgdHZ20tzcTGFh4ZcyvizLHDx4kNmzZ5OUlDSseERzz9AqdxuGBD7v9JASZ+byNbPJSI4Le+4eq5vGXjuZyTFkpYTvA+3s7KSlpWUY8UcC//jHP/jlL3/J7t27I9Iicc0113DmmWeyYcMG3G43drv9S+lBPQERJb4ovjxo4bbGxsZhVaitra2oqsqMGTOQJInDhw/z85//nPz8fHJzc4fpZ0bq+jweD1arlcrKShITExFFcVgINXi3GIm2ka6uLg4fPsyKFSu+FNILVT0ZDsF6qNr349FDDQWtt3XFihVTuvgJhYqKCq677jreeuutKXP28Mfg4CAFBQXU19cfdVPqKKLEF8VXGLIs88ADD/Dqq69y9dVX09LSEqBuk5SUpOcUNcm3cOo2E4XX6+XgwYPMmTNnmELHeEOoYy2bD0Z3dzcNDQ1H5YEfCpqR77Rp0yblaaidayx6qMGfm8PhoLy8nBUrVui53Uihra2NSy65hOeee46lS5eO/gsTQFFRETfccAOLFy+muLiYVatW8fDDD38pbTknIKLEF8VXFw6Hg02bNvHb3/522MNOVVX6+/tDqtsMDAwEqNtoYdR58+Yxc+bMMYfIPB4PRUVFIUlvNPg3Wfu3ZwSHUIN3i8G7ud7eXurq6r400lNVlfLycj1PG+mx/PVQNXK0Wq3YbDbi4uJ0STj/hcREtW1DYWBggLVr1/Lb3/6Wc845Z0rOGQqfffYZp556Kh999BGnnHIKt9xyC0lJSdx7770RGzMKHVHiGwteeuklNm/eTGVlJfv372f16tX6a/fffz9btmzBYDDwpz/9KaQJZUNDA+vWraOvr4+VK1fy17/+dcoT5VF8AX91m+AQakdHB0CAuo22W/RXt+nu7ubvf/87Z5111jBro6m4Pq1sPni36B9CBd+DeMGCBSQlJU3pA36s11lVVYUkSeTn5x+1cf3hcrk4ePAgixcvJj4+fsQdtr8eqn+R0lgXOi6Xi0svvZTrr7+edevWRfJt0dHRwamnnkpjYyMAH374If/1X//Fm2++GdFxowCixDc2VFZWIooiN954Iw899JBOfBUVFVxxxRW6rdF5551HdXX1sIl22WWX8Z3vfId169bxgx/8gIKCAl1tJYqjC03dprm5Oay6TWJiIq2trZx33nmcffbZk1K3mQhkWaazs5P6+noyMzP1nFlwCDU4HDjVFY41NTXIsszChQu/lDyUx+Ph4MGD5Ofnj+q2EayHGrzDDm5M9zcgFgQBRVG4/vrrWbVqFbfffvtReb9nnnkmTzzxBAsXLmTz5s3YbDYefPDBiI8bRZT4xoWzzz47gPjuv/9+AO666y4ALrjgAjZv3hxQkqyqKunp6XR0dCBJEnv37mXz5s3s2rXr6L+BKEZFb28vF154Id/5zneYNWvWMHWbmJiYAHWbefPmhVS3mQw0e6fCwsJhThehQqj+FZXhcmTjDZM2NDRgs9lYsmTJl0J6mkhAdnb2lLgfBDema5/d3/72N95++2194XD99dczb948/e8aSRQVFekVnXl5eTz11FOkpqZGdMwogKhyy+TQ2toaIFSbnZ1Na2trwDG9vb2kpKTouZtQx0Tx1UF5eTl333033/zmN4e9FkrdZvv27bq6jSiKZGdnD2vinzNnDpIkjYlABgYGwpIeEFAVGbwL8g+hao38PT09w0KowaQYHEJtbm5mcHCQZcuWfSmkp6oqZWVlzJgxY8osf4Ib0zWsWrWKOXPm8MEHH3DTTTfR2NjI7t276ejo4K9//euUjB0OhYWFUXPZrxhOOOI777zz9PyPP+67776wgrShdsXBD4qxHDMSLr/8cqqqqgDfTiAlJYWioqJhx0XdmqcGI6nfj0Xdxr814/3332fLli0B6jYaKfoTY3x8PIIgsG/fPnp6ejj33HMn5GkoCAImkwmTyRSy70yrQtW+Ojs7cTgcARJm2o4yPz8fh8MRkRDqSFBVlcrKShISEkJaPE013njjDV599VV27typ51WnClVVVVx++eX6z/X19dxzzz3DBCCi+OrghCO+d955Z9y/k52dTXNzs/5zS0vLMO3G6dOnY7FY8Hq9SJIU8piRsG3bNv3722+/fcRG2vfffz9qWvklQRAEjEYj8+fPZ/78+cNeD1a3KS8v5/XXX6ehoQG73Y4gCPT09HDZZZcxODg4YXWbkTCSFqqqqrS0tNDS0sLs2bMZGBigo6NjWFP6ZEOoo6Gurg5RFCMeZgTYv38/DzzwALt3755y0gNYuHChvkiVZZmsrCxdfzOKryZOOOKbCNauXcv69eu57bbbaGtro6amhpNPPjngGEEQOOecc9i+fTvr1q3j6aefnpCliaqqvPjii7z33ntTdflRHEWIoqg3fwfvKisqKrj88st5+OGHkWWZuro69u7dO0zdJtg1IysrK0DdZjLo6+ujvb2d1atXhySzYDupnp4evSk9nPnweKtQm5qasNvtRyXEWlNTw3/8x3+wY8cO0tLSIjoWwLvvvsu8efOYO3duxMeKYuKIFrf44ZVXXuFHP/oR3d3dpKSkUFhYqBen3HfffTz55JNIksQf//hHLrzwQgAuuuginnjiCWbNmkV9fb3ezrBixQqeffbZcTfh/uMf/+C2224LG8KMujUfu3j55ZdZtGgRixcvHvbaaOo2iqIwc+bMAB1UjSDHqm5jsVioqqpixYoVE2qzCQ6hagTpH0IN5bPov5Pt6OjQNUgjHVrt6uri4osv5oknnmDVqlURHUvDtddey8qVK9m4ceNRGS+KERGt6vwqYCz5xZtuuon8/Hxuv/32kOeYKrfmzZs385e//IX09HQAfvOb33DRRRcNO27nzp3ccsstyLLMhg0b2LRp07jHimLyUBSFtra2gEb+2trakOo2/nqo2iLp448/pru7mwsuuGBCecXRENxm4K/UorUZCIKA3W7XW0YiEULVYLPZ+Na3vsUvfvELfZEaabjdbmbNmkV5efmUFetEMSlEie9YgNfrJSsriwMHDowp4T8Zt+ax/K4syyxYsIA9e/aQnZ3NmjVreOGFF0LuWKL48jCauo0sy/T39/Ptb3+bpUuXTkjdZrLo7e3l0KFDzJ07NyCc6q/rGc5Id7zweDysX7+eSy65hGuvvTYC7yY0duzYwaOPPsru3buP2phRjIhoO8OxgHfeeYeTTjopLOkFuzXv3r2bu+++O2LXs3//fvLz8/UChHXr1rFjx44o8X3FIAgC06ZN4+STTx6We66vr+fb3/42//3f/40gCNTW1vL6669TX19Pe3s7giCQmZk5bLeYk5MzZSLcNpuN6upqVq5cGbK4RJblAOeMrq6uYSHUUD6LoUKliqJw2223ccopp/D9739/0tc+HrzwwgtcccUVR3XMKCaGKPF9hbB169ZhE6etrY0NGzbw1ltvTblb8yOPPMIzzzzD6tWr+d3vfjesqba1tTVArDg7O5t9+/ZNeLwojj5qa2t56qmnQua4Qqnb7N27l+eee47Dhw/jdruZNm1aQKHNeNVtnE4npaWlLF26NGxFpcFgID4+PqRwc3AIdWBggPb2dhwOB6qq6kotb7/9NikpKVRWVuJ2u/n5z39+VHsT7XY7e/bs4fHHHz9qY0YxcURDnccxRsopnnrqqUyfPh1BEPjFL35Be3s7Tz75ZMBx/7+9ew9p8o3iAP6dLH9slP4hKk4pb2FLRV0XrQwp2DIjTS3tKtEoDCQlCCFBKEgLukkXtSSzUS00hxlbt5EUXhhJM7QLswtaqKuF6RJ05vv7IxwtZ2l759SdDwjb3gfPI4LH93nec57Kykrcv38fZWVlAACZTAaNRoNz585Nei6HDh1CbW0tXF1dERQUhPLycqtnklGd4vQxMjICg8FgdQnVaDSO6W4TGBiIoKAgc3ebz58/Q6FQYMuWLXbrVDJayF9TUwO1Wo2WlhYIBAJzA/O0tDQcPHjQLrHJtEdLnc5oojWLe/futdrBZCL1ixMlFotRWFgILpeL3NxcFBYW4sSJE1bHUp3i9ODi4gJPT094enqOOTHcWneb27dv4+3bt9Dr9eBwOOjt7cXKlSvBMIy5QfhkuttMxGinFh8fH3z8+BHNzc3m+sXBwUEYjUZW4pDZhRKfk+rq6jIfNKpQKKyeR7Zs2TLodDq8f/8evr6+kMvluHHjxj/Fk0gk5tcxMTGoqqr6t4mTaeFP3W1MJpO5WXtUVBTa29tRV1eHK1euoLOzE8PDw3/tbjMZWq0W+fn5uHfvnkXR/n///cfKmX5nzpxBWVkZOBwOwsPDUV5ebpcnY8nUoaVOJ7Vr1y5otVpwOBz4+/ujtLQUPj4+FnuKAKBUKpGTk4MfP35gz549yMvLszn2xo0bkZ6ejp07d465RnWKM19XVxdqa2vH/d393t1mdPl0tLsNn883J0N/f39zI2lr3W06OjqQlpaGyspKhISEsP6zfPr0CbGxsXj58iV4PB7S0tKQkJCA3bt3sx6LsI7KGYj9TaRO8dixY3j27Bmqq6ut/mfPRp3i32oPBwcHkZGRgebmZnh4eODWrVvw9/efVAxiHwzDoK+vb0wR/2h3Gy6Xa+5u4+3tjYqKCly+fBmrVq2yy3xGG9S3tLTAzc0NmzZtwoEDByxWMMi0RYmPOF5FRQVKSkqgVqvB5/P/Ov5f6hQnUnt48eJFvHjxAiUlJZDL5VAoFBa9Usn09Ht3m8bGRvB4PFZWIf6kqKgIeXl54PF4kEgkuH79ul3jEdaMv2bOMMyfvghhhUqlYoRCIaPX68cdYzQamb6+PvPrFStWMCqValJxGhoaGIlEYn5fUFDAFBQUWIyRSCRMQ0MDwzAMYzKZGA8PD2ZkZGRScYhz+Pr1K7NmzRpGr9czQ0NDTFJSEiOTyRw9LTIx4+a2qTuHhDi1rKws9Pf3QywWIzIyEpmZmQB+Lm2Otkrr6elBbGwsIiIisHz5cmzYsGHSdYrWag9/Pxfx1zFcLhfu7u4wGAy2/Hhklnr06BECAgLg6emJOXPmICUlBQ0NDY6eFrERPdVJpkR7e7vVzwUCgflBmsDAQLS0tNgUh5mCsxN/19nZiYyMDHR3d8PFxQX79u1Ddna2xZi6ujokJSUhICAAAJCSkmLXrjuEHfPnz0dTUxMGBgbA4/GgVquxdOlSR0+L2IgSH5lVJlJ7ODrGz88Pw8PD+Pbt25hTzieDy+Xi1KlTEIlE6O/vx5IlSyAWi8e0dlu9ejXu3r37z3HI1IuOjsbmzZshEonA5XIRFRVFTxrPArTUSWaVX+dQpEoAAAPpSURBVGsPh4aGIJfLkZiYaDEmMTERFRUVAICqqiqsXbvWpjs+Hx8fiEQiAMC8efMgFArHLK+SmevIkSN4/fo1WltbIZPJWKkNJI5FiY/MKlwuF+fPn8e6desgFAqRlpaG0NBQ5Ofn486dOwAAqVQKg8GA4OBgnD59GsePH2ct/ocPH/D8+XNER0ePudbY2IiIiAisX78ebW1trMUktikqKkJYWBhCQ0Nx9uxZR0+HTAEqZyCEJUajEXFxccjLy0NKSorFtb6+Pri4uGDu3LlQKpXIzs6GTqezKd7f+poyDIPs7GwolUrw+XxcvXrVfGdKfmptbcXWrVuh0Wjg6uqK+Ph4FBcXY+HChY6eGrHduMs4dMdHCAtMJhNSU1OxY8eOMUkPANzc3MzttBISEmAymfDlyxeb4z5+/BhardZqM2+VSgWdTgedTodLly5h//79NsebbV69eoWYmBjw+XxwuVzExcVBoVA4elrEzijxEWIjhmEglUohFArHPQmgu7vb/DSpRqPByMgIPDw87DqvmpoaZGRkgMPhICYmBr29vejq6rJrzJkmLCwMT548gcFgwMDAAJRKpcXDUWR2oqc6CbFRfX09ZDIZwsPDERkZCQAoKChAR0cHACAzMxNVVVUoLi42nx8nl8ttPqGAw+FAIpGM29d0vJrG0ebkBBAKhcjNzYVYLDY33eZy6c/ibEe/YUJsFBsba7U28FdZWVnIyspiNW59fb1FX9NFixZZ9DVlu17xzZs3SE9PN79/9+4djh49ipycHPNnM7FeUSqVQiqVAgAOHz4MPz8/B8+I2BslPkJmqNH6RC8vLyQnJ0Oj0VgkPjbPUwSAkJAQaLVaAD97ovr6+iI5OXnMuJlWr6jX6+Hl5YWOjg5UV1ejsbHR0VMidkZ7fITMQN+/f0d/f7/59YMHD8acqZiYmIhr166BYRg0NTXB3d2dtWVOtVqNoKAgLFiwgJXv50ipqalYvHgxNm7ciAsXLtjttHgyfdAdHyEzUE9Pj/lua3h4GNu3b0d8fDxKSkoA/NxXTEhIgFKpRHBwMPh8PsrLy1mLL5fLsW3bNqvXRusVBQIBTp48idDQUNbi2sPTp08dPQUyxaiOjxAyKUNDQxAIBGhra4O3t7fFNXvUKxLyj6iOjxDCDpVKBZFINCbpAfarVySETZT4CCGTcvPmzXGXOR1Rr0jIZNEeHyFkwgYGBvDw4UOUlpaaP/t1X9Ee9YqEsI32+AghhMxGtMdHCCGEAJT4CCGEOJm/7fHR4jwhhJBZhe74CCGEOBVKfIQQQpwKJT5CCCFOhRIfIYQQp0KJjxBCiFOhxEcIIcSp/A+JEvuZK2slywAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# 绘制抽样样本的二维正态分布\n",
    "from mpl_toolkits.mplot3d import Axes3D\n",
    "fig = plt.figure()\n",
    "ax = Axes3D(fig, rect=[0, 0, 1, 1], elev=30, azim=20)\n",
    "ax.scatter(x_res, y_res, z_res,marker='o')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.7"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
